Dynamic fundus image registration method and system based on space-time coupling

Abstract

The application relates to the technical field of image processing, and discloses a dynamic fundus image registration method and system based on space-time coupling, which comprises the following steps: positioning the light source brightness energy center coordinates of a current processing frame of a laser speckle fundus image sequence and calculating the time sliding window gray scale gradient fluctuation rate; performing coupling multiplication operation on a space confidence parameter and a time confidence parameter to generate a pixel-level dynamic confidence tensor; calculating a global motion divergence parameter and comparing the global motion divergence parameter with a demarcation threshold value to determine a motion state interval of the current processing frame; triggering a corresponding registration processing path according to the motion state interval; finally generating an updated dynamic synthesis reference frame; and reversely adjusting the demarcation threshold value based on the updated dynamic synthesis reference frame; the application realizes pixel-level isolation of an interference area and closed-loop feedback of a registration strategy, and improves the robustness of dynamic fundus image registration.

Description

Technical Field

[0001] This invention relates to the field of image processing technology, and more specifically, to a method and system for dynamic fundus image registration based on spatiotemporal coupling. Background Technology

[0002] In the field of fundus medical image analysis, laser speckle fundus imaging technology is a key means of observing retinal hemodynamics and microvascular morphology. In order to obtain blood flow perfusion maps with high signal-to-noise ratio and eliminate speckle noise in single-frame images, it is usually necessary to perform multi-frame registration and weighted fusion operations on continuously acquired fundus image sequences to achieve precise alignment of retinal vascular topology in different frames in the same spatial coordinate system.

[0003] In the prior art, patent CN101937565B discloses a dynamic image registration method based on the trajectory of a moving target. This method extracts the connected regions of the moving target and generates a trajectory by associating it with previous and subsequent frames. It uses chain code to describe the trajectory and performs matching to obtain temporal and spatial registration parameters, thus achieving inter-frame registration under conditions of fixed or changing backgrounds. Patent application CN116245923A discloses an image registration method and apparatus, as well as an image processing method. It utilizes a deep learning model containing deformable convolutional layers and fully connected layers to determine deformation data by capturing the feature point matching relationship between the fundus image to be corrected and the reference image. This aims to reduce the impact of resolution differences, geometric distortions, and interference such as occlusion and blurring on the registration accuracy.

[0004] However, in actual clinical scenarios of laser speckle fundus imaging, when the subject has excessively long eyelashes or intermittent eyelid micro-closure, gray stripes that flicker frequently with breathing or eye movements will be introduced at the image edges. Because these interfering elements are highly similar to peripheral blood vessels in their gray-scale gradient characteristics, and because the radial distortion caused by the light source illuminating the bowl-shaped structure of the eyeball is unevenly distributed, traditional global registration algorithms or fixed feature point matching mechanisms are prone to misidentifying eyelash edges or skin reflection spots as vascular anchor points when processing such sequences. This misidentification leads to incorrect deformation compensation during registration, resulting in vascular texture ghosting or "phantom" phenomena in the dynamically synthesized reference image. If the system fails to identify and isolate these low-quality frames caused by physiological intermittent eye movements in real time, erroneous information will continuously pollute the reference database, ultimately causing directional deviations in the velocity vector field calculations during subsequent quantitative hemodynamic analysis, failing to accurately reflect the true state of retinal circulation. Summary of the Invention

[0005] The core of the aforementioned problem lies in the nonlinear evolution of "spatiotemporal coupling interference" during fundus image acquisition. Laser speckle imaging is highly sensitive to phase changes; even minute displacements of eyelashes at the edge of the field of view can cause dramatic fluctuations in the grayscale gradient of target pixels over time, exhibiting statistically high coefficients of variation. Simultaneously, the posterior pole of the eyeball, due to its bowl-shaped depression, naturally exhibits non-rigid geometric distortions in the spatial dimension in areas far from the light source's brightness energy center. When these two types of interference overlap, a single pixel grayscale value or a static feature extraction operator cannot distinguish between "stable vascular texture" and "fluctuating eyelash projection." If the system only employs a uniform registration strategy, high-quality information in a sub-pixel tremor state will overlap with invalid information in a blink-locked state during the fusion process. This lack of discriminative processing logic prevents the effective suppression of global motion divergence peaks generated during large eye movements, resulting in the loss of constraint on spatial center deviation distance and temporal gradient fluctuation rate at the pixel level. Because there is no feedback loop between the registration execution stage and the state determination stage, the image information entropy of the reference frame will decay after the interference frame enters. If this decay cannot act in reverse to adjust the motion state boundary threshold, it will cause the entire registration closed loop to collapse, making it impossible for the image reference library to maintain sub-pixel level registration accuracy under complex physiological environments.

[0006] To overcome the aforementioned shortcomings of existing technologies, this invention provides a dynamic fundus image registration method and system based on spatiotemporal coupling. By coupling the spatial constraints driven by the energy center of the coupled light source with the temporal constraints driven by grayscale gradient fluctuations, a pixel-level dynamic confidence tensor is constructed to achieve flexible isolation of interference regions. The system divides motion state intervals based on global motion divergence parameters and triggers differentiated processing paths. It uses an updated dynamically synthesized reference frame to inversely adjust the judgment threshold, maintaining the robustness of the registration process and the clarity of the fused image under complex eye-tracking conditions.

[0007] To achieve the above objectives, the present invention provides the following technical solution:

[0008] A dynamic fundus image registration method based on spatiotemporal coupling includes:

[0009] For the current processing frame in the laser speckle fundus image sequence, the coordinates of the light source brightness and energy center are located and the gray-level gradient fluctuation rate of the time sliding window is calculated. The spatial confidence parameter driven by the light source brightness and energy center coordinates and the temporal confidence parameter driven by the gray-level gradient fluctuation rate are coupled and multiplied to generate the pixel-level dynamic confidence tensor of the current processing frame. The global motion divergence parameter is calculated and compared with the first state boundary threshold and the second state boundary threshold to determine the motion state interval to which the current processing frame belongs.

[0010] The first frame of the laser speckle fundus image sequence is used to initialize the dynamic synthesis reference frame. The corresponding registration processing path is triggered according to the motion state interval to which the current processing frame belongs. This includes performing fine topological registration and outputting the registered current frame when the current processing frame is in the first state interval. The fusion ratio is calculated based on the pixel-level dynamic confidence tensor. The registered current frame and the dynamic synthesis reference frame are weighted and fused according to the fusion ratio to obtain the updated dynamic synthesis reference frame. The first state boundary threshold and the second state boundary threshold are adjusted based on the updated dynamic synthesis reference frame.

[0011] The method for locating the center coordinates of the light source's brightness energy includes:

[0012] Read the gray values ​​of all pixels in the current processing frame of the laser speckle fundus image sequence, and perform weighted centroid calculation on the gray values ​​of all pixels in the current processing frame to obtain the coordinates of the light source brightness energy center of the current processing frame.

[0013] The method for calculating the time-sliding window gray-scale gradient volatility includes:

[0014] Using the current processing frame as the last frame, the preceding consecutive frames are taken to form a time sliding window. For each target pixel within the time sliding window, the local gray-level gradient magnitude sequence is extracted. The standard deviation and mean of the local gray-level gradient magnitude sequence are calculated. The gray-level gradient fluctuation rate of the target pixel is obtained based on the ratio of the standard deviation to the mean of the gray-level gradient magnitude sequence.

[0015] The method for generating the pixel-level dynamic confidence tensor of the current processing frame includes:

[0016] The spatial confidence parameters of the target pixel are constructed based on the coordinates of the light source brightness energy center, and the temporal confidence parameters of the target pixel are constructed based on the gray-scale gradient fluctuation rate.

[0017] The spatial confidence parameter and temporal confidence parameter of the target pixel are coupled and multiplied to obtain the pixel-level dynamic confidence value of the target pixel. The pixel-level dynamic confidence values ​​of all target pixels in the current processing frame are then combined to form the pixel-level dynamic confidence tensor of the current processing frame.

[0018] The method for constructing the spatial confidence parameters of the target pixel based on the coordinates of the light source brightness energy center includes:

[0019] The center offset distance of the target pixel is obtained by calculating the Euclidean distance between the target pixel and the center coordinates of the light source brightness energy in the current processing frame. The spatial confidence parameter of the target pixel is constructed based on the center offset distance.

[0020] The calculation method for the global motion divergence parameter includes:

[0021] The current processing frame is divided into multiple non-overlapping rectangular sub-blocks. For each rectangular sub-block, the matching position with the smallest grayscale difference in the previous frame is searched to obtain the local displacement vector. The arithmetic mean of the magnitudes of all local displacement vectors is calculated to obtain the displacement magnitude mean. The standard deviation of the magnitudes of all local displacement vectors is calculated to obtain the displacement magnitude standard deviation. The displacement magnitude mean and the displacement magnitude standard deviation are added together to obtain the global motion divergence parameter.

[0022] The method for determining the motion state interval to which the current processing frame belongs includes:

[0023] When the value of the global motion divergence parameter is less than the first state boundary threshold, the current processing frame is determined to be in the first state interval, which represents the subpixel micro-tremor state; when the value of the global motion divergence parameter is greater than or equal to the first state boundary threshold and less than the second state boundary threshold, the current processing frame is determined to be in the second state interval, which represents the large-amplitude eye movement transition state; when the value of the global motion divergence parameter is greater than or equal to the second state boundary threshold, the current processing frame is determined to be in the third state interval, which represents the blink lock state; wherein, the first state boundary threshold is less than the second state boundary threshold.

[0024] The method for performing fine-grained topology registration includes:

[0025] Pixels with pixel-level dynamic confidence values ​​greater than a set confidence threshold are extracted from the pixel-level dynamic confidence tensor of the current processing frame to form a high-confidence pixel set; the high-confidence pixel set is used to fit the non-rigid deformation field between the current processing frame and the dynamically synthesized reference frame, and the current processing frame is resampled according to the non-rigid deformation field to output the registered current frame.

[0026] The step of triggering the corresponding registration processing path based on the motion state interval to which the current processing frame belongs also includes:

[0027] When the motion state interval determination result indicates that the current processing frame is in the second state interval, the rigid anchor suspension mode is triggered; when the motion state interval determination result indicates that the current processing frame is in the third state interval, the trajectory deduction skip mode is triggered.

[0028] The rigid anchor point suspension mode includes extracting the main vessel intersection anchor points from the pixel-level dynamic confidence tensor of the current processing frame, calculating the median of the coordinate offset between the main vessel intersection anchor points and the dynamic synthesis reference frame, outputting the whole frame coarse rigid offset vector of the current processing frame, recording the whole frame coarse rigid offset vector as the whole frame unified rigid offset vector of the current processing frame to the motion trajectory history sequence, and prohibiting the current processing frame from participating in pixel-level fusion.

[0029] A spatiotemporal coupling-based dynamic fundus image registration system is used to implement the above-mentioned spatiotemporal coupling-based dynamic fundus image registration method. The system includes:

[0030] Confidence Tensor Generation Module: Performs source brightness and energy center coordinate localization and time sliding window grayscale gradient fluctuation rate calculation on the current processing frame in the laser speckle fundus image sequence. Couples and multiplies the spatial confidence parameter driven by the source brightness and energy center coordinate with the time confidence parameter driven by the grayscale gradient fluctuation rate to generate the pixel-level dynamic confidence tensor of the current processing frame.

[0031] Motion state determination module: Performs block matching calculation on the current processing frame and the previous frame to obtain global motion divergence parameters, and compares the global motion divergence parameters with the first state boundary threshold and the second state boundary threshold to determine the motion state interval to which the current processing frame belongs;

[0032] Reference frame registration trigger module: Initializes the dynamic synthesis reference frame with the first frame of the laser speckle fundus image sequence, and triggers the corresponding registration processing path according to the motion state interval to which the current processing frame belongs;

[0033] Fine registration and fusion module: When the current processing frame is in the first state interval, fine topological registration is performed and the registered current frame is output. The fusion ratio is calculated based on the pixel-level dynamic confidence tensor. The registered current frame and the dynamic synthesis reference frame are weighted and fused according to the fusion ratio to obtain the updated dynamic synthesis reference frame. The first state boundary threshold and the second state boundary threshold are adjusted based on the updated dynamic synthesis reference frame.

[0034] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0035] This invention addresses the feature extraction mismatch problem caused by local fundus distortion, eyelash occlusion, and intermittent eye movements during fundus image acquisition. It constructs a pixel-level dynamic confidence tensor by coupling and multiplying spatial constraints driven by the light source brightness energy center with temporal constraints driven by grayscale gradient volatility. This achieves adaptive and flexible isolation of static field-of-view occlusion and dynamic abrupt interference at the pixel level. Simultaneously, a global motion divergence parameter is introduced and combined with a dynamic boundary threshold to divide motion state intervals. Differentiated registration processing paths are triggered based on the intensity of eye movements, avoiding systemic bias caused by forced registration of low-quality interference frames. The pixel-level dynamic confidence tensor is used as the calculation benchmark for the fusion ratio to weighted update the dynamically synthesized benchmark frame. The updated benchmark frame is then used to inversely adjust the boundary threshold of the motion state. A system-level dynamic closed-loop feedback network is established between registration execution and parameter determination, effectively preventing the continuous contamination of the image benchmark library by intermittent physiological interference and improving the registration accuracy and fusion reliability of temporal images under complex eye movement conditions. Attached Figure Description

[0036] 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.

[0037] Figure 1 This is a flowchart of the dynamic fundus image registration method based on spatiotemporal coupling in this invention;

[0038] Figure 2 This is a schematic diagram of the gray-level gradient amplitude sequence within the time sliding window in this invention;

[0039] Figure 3 This is a schematic diagram illustrating the division of the local displacement vector motion state interval in the block matching method of this invention.

[0040] Figure 4 This is a flowchart of the conditional sub-control criteria driven by the motion state interval in this invention.

[0041] Figure 5 This is a schematic diagram of the extraction of the main blood vessel confluence anchor points in this invention;

[0042] Figure 6 This is a functional block diagram of the dynamic fundus image registration system based on spatiotemporal coupling in this invention. Detailed Implementation

[0043] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0044] Example 1:

[0045] Please see Figure 1 As shown, this embodiment provides a dynamic fundus image registration method based on spatiotemporal coupling, including:

[0046] Step S10: Perform source brightness and energy center coordinate localization and time sliding window grayscale gradient fluctuation rate calculation on the current processing frame in the laser speckle fundus image sequence. Couple the spatial confidence parameter driven by the source brightness and energy center coordinate with the time confidence parameter driven by the grayscale gradient fluctuation rate to generate the pixel-level dynamic confidence tensor of the current processing frame. Calculate the global motion divergence parameter. Compare the global motion divergence parameter with the first state boundary threshold and the second state boundary threshold to determine the motion state interval to which the current processing frame belongs.

[0047] During clinical acquisition, laser speckle fundus image sequences exhibit persistent moving gray stripes formed by eyelash projections and large white, informationless areas formed by saturated reflections of eyelid skin in the image edge regions. These interfering elements resemble real blood vessel features in terms of gray-scale gradient morphology. Once extracted as feature points, they will mismatch with blood vessel feature points in the reference frame. Simultaneously, fundus image sequences mostly exhibit only sub-pixel-level physiological micro-tremors between adjacent frames, but intermittently insert low-quality frames caused by blinking or large eye movements. If all frames are indiscriminately sent to a unified registration process, the mismatch results of low-quality frames will be passed on to subsequent multi-frame fusion stages, introducing systematic bias. Step S10 constructs a pixel-level dynamic confidence tensor to suppress the weights of interfering regions and calculates global motion divergence parameters and divides motion state intervals to achieve hierarchical management of frame quality. This provides two types of inputs for multimodal conditional subdivision and registration in step S20: pixel-wise reliability metrics and frame-wise motion state labels.

[0048] Further, step S10 includes:

[0049] Step S11: Obtain the laser speckle fundus image sequence; perform weighted centroid calculation on the gray values ​​of all pixels in the current processing frame of the laser speckle fundus image sequence to obtain the coordinates of the light source brightness energy center of the current processing frame; calculate the Euclidean distance between the target pixel and the coordinates of the light source brightness energy center in the current processing frame to obtain the center deviation distance of the target pixel; and construct the spatial confidence parameter of the target pixel based on the center deviation distance.

[0050] After acquiring the laser speckle fundus image sequence, the grayscale matrix of the current processing frame is read from the sequence. This grayscale matrix consists of the grayscale values ​​of all pixels in the current processing frame. The physical meaning of the light source brightness energy center coordinates is the weighted centroid position of the grayscale value distribution in the current processing frame, representing the projection center of the laser source on the fundus imaging plane. To calculate the light source brightness energy center coordinates, the row coordinate of each pixel in the current processing frame is multiplied by its corresponding grayscale value. The sum of these products for all pixels is then divided by the total grayscale value of all pixels in the current processing frame. The quotient is the row component of the light source brightness energy center coordinates. Similarly, the column coordinate of each pixel in the current processing frame is multiplied by its corresponding grayscale value. The sum of these products for all pixels is then divided by the total grayscale value of all pixels in the current processing frame. The quotient is the column component of the light source brightness energy center coordinates. The weighted centroid calculation described above belongs to the gray-scale centroid method in the field of image processing. It uses gray values ​​as weights for spatial location, so that pixels with higher brightness contribute more to the centroid location. Since the brightness produced in the central region of the imaging surface after a laser light source illuminates the fundus is much higher than that in the edge region, the gray-scale centroid method can locate the coordinates of the light source's brightness energy center to the core region illuminated by the light source, rather than being offset by interference from dark areas at the edges or areas obscured by eyelashes.

[0051] Based on the coordinates of the light source's brightness energy center, for any target pixel in the current processing frame, the square root of the square difference between the row coordinates of the target pixel and the row component of the light source's brightness energy center coordinates, plus the square difference between the column coordinates of the target pixel and the column component of the light source's brightness energy center coordinates, is taken. The resulting value is the center offset distance of the target pixel, expressed in pixels. A spatial confidence parameter for the target pixel is constructed based on this center offset distance. The formula for calculating the spatial confidence parameter is as follows:

[0052] ;

[0053] In the formula, This represents the spatial confidence parameter of the target pixel. Indicates the first attenuation base. This indicates the center offset distance of the target pixel. This represents the first normalization factor. The first decay base. A constant greater than 1, used to control the rate basis for the decay of the spatial confidence parameter as the center deviation distance increases. First normalization factor. It is equal to the reciprocal of the diagonal pixel length of the current processing frame. The diagonal pixel length of the current processing frame is the square root of the sum of the squares of the number of pixels in the row direction and the squares of the number of pixels in the column direction of the current processing frame. The first normalization factor is set to the reciprocal of the diagonal pixel length of the current processing frame, so that the center deviation distance in the exponent is normalized to the range of 0 to 1. When the target pixel is located at the center coordinate of the light source brightness energy, the center deviation distance is 0 and the spatial confidence parameter is 1. When the target pixel is located at the diagonal distance of the current processing frame, the normalized exponent term is -1 and the spatial confidence parameter is the -1 power of the first attenuation base. The determination of the first attenuation base needs to consider two factors: the distortion distribution characteristics caused by the eye's bowl-shaped structure on the imaging surface and the radial attenuation law of the light source illumination. If the first attenuation base is too small, only slightly greater than 1, the spatial confidence parameter decays too slowly as the center deviation distance increases. The edge region far from the center coordinates of the light source brightness energy still maintains a high spatial confidence parameter value, resulting in the eyelash interference pixels in the edge region not being fully suppressed and still being included in the high-confidence pixel set for matching in the subsequent step S21. If the first attenuation base is too large, the spatial confidence parameter decays to an extremely low value after a very short distance from the center coordinates of the light source brightness energy, causing the pixels in the optic disc and surrounding blood vessel regions to be unable to participate in feature matching due to the low spatial confidence parameter, while the optic disc region is precisely the core region of interest for registration. In actual determination, the ratio of the effective coverage radius of the light source on the fundus imaging surface to the diagonal pixel length of the current processing frame is used as a reference, and a value that causes the spatial confidence parameter to decay to the range of 0.3 to 0.5 at the effective coverage radius is selected as the first attenuation base. The effective coverage radius refers to the radial distance from the peak value at the center coordinates of the light source's brightness energy to half the peak value after the laser light source illuminates the fundus imaging surface. This distance reflects the actual useful boundary of the illuminated area. In a typical laser speckle fundus imaging system, the effective coverage radius usually accounts for 30% to 45% of the diagonal pixel length of the current processing frame, meaning that the area effectively illuminated by the light source occupies about one-third to nearly half of the entire image area.

[0054] The spatial confidence parameter is constructed using an exponential decay method instead of setting rigid region boundaries. This allows the spatial confidence parameter to decrease continuously and smoothly outward from the center coordinates of the light source's brightness energy, without any abrupt cutoff points. Pixels in the edge transition region receive intermediate weights between 0 and 1, rather than being crudely retained or discarded. This continuous weight allocation mechanism is consistent with the physical law of continuously increasing distortion caused by the bowl-shaped structure of the eyeball. This allows the spatial confidence parameter to implicitly quantify the degree of image distortion solely through the topological distance between the pixel and the center coordinates of the light source's brightness energy, without relying on any edge detection or region segmentation algorithms. The spatial confidence parameter output in step S11 will be coupled with the temporal confidence parameter in step S13. Without the spatial dimension constraint provided in step S11, the pixel-level dynamic confidence tensor generated in subsequent step S13 will only have temporal dimension discrimination ability and will be unable to identify pixels that, although their gray-level gradients are stable in the temporal dimension, are actually located in high-distortion edge regions of the image. If these pixels participate in the non-rigid deformation field fitting in step S21, they will introduce systematic deformation errors.

[0055] Step S12: Using the current processing frame as the last frame, take the preceding consecutive frames to form a time sliding window. For each target pixel within the time sliding window, extract the local gray-level gradient magnitude sequence, calculate the standard deviation and mean of the local gray-level gradient magnitude sequence, obtain the gray-level gradient volatility of the target pixel based on the ratio of the standard deviation to the mean of the gray-level gradient magnitude sequence, and construct the time confidence parameter of the target pixel based on the gray-level gradient volatility.

[0056] See Figure 2 This is a schematic diagram of the gray-level gradient magnitude sequence within a time sliding window provided in an embodiment of this application, wherein... Figure 2 (a) illustrates a time-sliding window that takes the current processing frame as the last frame and includes preceding consecutive historical frames. Figure 2(b) illustrates the sequence curve of the grayscale gradient magnitude corresponding to a target pixel within the window as it changes over time, showing the abrupt change in grayscale gradient magnitude when the eyelash projection sweeps across the pixel position. In step S12, taking the current processing frame as the last frame, the first consecutive number of historical frames are taken backwards, forming a time sliding window together with the current processing frame. Here, the target pixel refers to a pixel position in the current processing frame, and the pixel at the same row and column coordinate position is taken in each frame within the time sliding window. Determining the first number requires considering two factors: the frame rate of the laser speckle fundus image sequence and the periodic characteristics of physiological micro-tremors of the eye. If the first number is too small, the time span covered by the time sliding window is insufficient to encompass the complete process of an eyelash projection sweeping across a pixel location, resulting in the gray-level gradient fluctuation rate failing to adequately reflect the temporal instability of interfering pixels. This makes it impossible for interfering pixels and vascular pixels to produce sufficient distinction in terms of temporal confidence parameters. If the first number is too large, the time span covered by the time sliding window exceeds the period of eye micro-tremors, and the cumulative displacement between frames within the window exceeds the sub-pixel level. The gray-level gradient changes of vascular pixels caused by displacement are incorrectly included in the fluctuation rate, causing the temporal confidence parameter of the real vascular pixels to be unreasonably suppressed. In practice, the frame rate of the laser speckle fundus image sequence multiplied by the typical single oscillation period of eye micro-tremors is used as a reference, and the number of frames that allow the time sliding window to cover exactly a certain number (e.g., 2-4) of micro-tremor oscillation periods is selected as the first number.

[0057] When the current processing frame is located at the beginning of the laser speckle fundus image sequence, and the number of consecutive historical frames available before the current processing frame is less than a first number, the total number of consecutive historical frames actually available before the current processing frame is defined as the actual number of usable frames. The first number is replaced by the actual number of usable frames to construct a time sliding window. This time sliding window consists of the actual number of usable historical frames and the current processing frame. When the actual number of usable frames is zero, i.e., the current processing frame is the first frame of the laser speckle fundus image sequence, the time sliding window contains only the current processing frame itself. The length of the grayscale gradient magnitude sequence is 1, the standard deviation is zero, the grayscale gradient volatility is zero, and the time confidence parameter is set to 1, indicating that no time-dimension confidence penalty is applied to the first frame when there is no historical frame reference. When the actual number of available frames is greater than zero but less than a first quantity, the length of the gray-level gradient magnitude sequence is equal to the actual number of available frames plus 1. Based on this length, the gray-level gradient volatility is calculated using the same ratio of standard deviation to mean as the full window. The resulting gray-level gradient volatility is still used to construct the temporal confidence parameter of the target pixel. Because the actual number of available frames is less than the first quantity, the time span covered by the time sliding window is shorter than expected. Insufficient statistical sample size for the gray-level gradient volatility may lead to statistical biases of higher or lower gray-level gradient volatility for individual pixels. However, as the laser speckle fundus image sequence progresses frame by frame, the number of consecutive historical frames that can be obtained before the current processing frame continues to increase until it reaches the first quantity. The time sliding window then returns to the full window length, and the statistical reliability of the gray-level gradient volatility returns to a normal level.

[0058] For each frame within a time-sliding window, the local gray-level gradient magnitude of each pixel in that frame is calculated. The gray-level gradient magnitude sequence is the result of chronologically arranging the local gray-level gradient magnitudes of the same target pixel in each frame within the window. The calculation of the local gray-level gradient magnitude uses the Sobel operator to perform convolution operations in both the horizontal and vertical directions. The square root of the sum of the squares of the horizontal and vertical convolution results yields the local gray-level gradient magnitude of the corresponding pixel in the corresponding frame. The Sobel operator is a differential operator used for edge detection in image processing. It estimates the rate of gray-level change in an image in two orthogonal directions by weighted differencing of the gray-level values ​​in the pixel's neighborhood. After performing the above operation on all frames in the time sliding window, the target pixel obtains a gray-level gradient magnitude sequence composed of local gray-level gradient magnitudes arranged in chronological order. When the number of consecutive historical frames that can be obtained before the current processing frame is greater than or equal to a first number, the length of the gray-level gradient magnitude sequence is equal to the first number plus 1; when the number of consecutive historical frames that can be obtained before the current processing frame is less than the first number, the length of the gray-level gradient magnitude sequence is equal to the actual number of consecutive historical frames that can be obtained before the current processing frame plus 1.

[0059] Based on the gray-level gradient amplitude sequence of the target pixel, the arithmetic mean of all elements in the gray-level gradient amplitude sequence and the standard deviation of all elements relative to the arithmetic mean are calculated. The gray-level gradient volatility of the target pixel is obtained by dividing the standard deviation by the arithmetic mean. Gray-level gradient volatility is the coefficient of variation, a dimensionless indicator in statistics that measures the degree of data dispersion. When the target pixel corresponds to the edge of a real blood vessel, the displacement of the blood vessel within the sub-pixel micro-twitch range has a very limited impact on the local gray-level gradient amplitude. The elements in the gray-level gradient amplitude sequence maintain similar values, with a small standard deviation and a large mean, resulting in a low gray-level gradient volatility. When the target pixel corresponds to the area swept by the eyelash projection, such as... Figure 2 (a) Within the time sliding window, frames where the eyelash projection covers the target pixel's location occur. The movement of the eyelash causes the target pixel to be occluded by the eyelash in some frames, resulting in a sudden increase in the grayscale gradient magnitude, while in other frames it is not occluded and the grayscale gradient magnitude decreases. Figure 2 (b) The gray-level gradient amplitude sequence shows a clear single peak, a large standard deviation, and a high gray-level gradient volatility.

[0060] Based on the grayscale gradient volatility of the target pixel, a temporal confidence parameter for the target pixel is constructed. The formula for calculating the temporal confidence parameter is as follows:

[0061] ;

[0062] In the formula, This represents the time confidence parameter of the target pixel. This represents the grayscale gradient fluctuation rate of the target pixel. Indicates the second normalization factor. This indicates taking the larger of the two values. Second normalization factor. A magnification factor greater than 1 is used to map the numerical range of gray-level gradient volatility to a scale that provides sufficient differentiation for the time confidence parameter between 0 and 1. The determination of the second normalization factor needs to consider the baseline gray-level gradient volatility level caused by the inherent speckle noise in the laser speckle image: if the second normalization factor is too small, only pixels with extremely high gray-level gradient volatility will cause the time confidence parameter to drop to a low value, resulting in some interfering pixels with moderate volatility still maintaining a high time confidence parameter and not being suppressed; if the second normalization factor is too large, normal baseline volatility caused by speckle noise can cause the time confidence parameter of large-area pixels to drop to zero, resulting in pixels in the effective blood vessel area being incorrectly suppressed. In actual determination, the mean of the gray-level gradient volatility of pixels in the central region without eyelash obstruction is selected from the calibration sequence as the baseline volatility reference, and a value that makes the time confidence parameter corresponding to the baseline volatility within the range of 0.7 to 0.9 is selected as the second normalization factor. The calibration sequence refers to a short sequence of fundus images acquired before formal acquisition, requiring the test subject to maintain a stable gaze. It typically contains 50 to 200 frames and is used to obtain basic parameter statistics for the current test subject and the current imaging system combination. Baseline volatility refers to the mean gray-level gradient volatility of pixels in the unobstructed central region of the calibration sequence within a time sliding window. This value reflects the normal gray-level gradient volatility level caused only by laser speckle noise and weak physiological micro-tremors, without any external interference factors. In typical laser speckle fundus imaging systems, the baseline volatility value is usually in the range of 0.05 to 0.15. Taking a larger value ensures that when the gray-level gradient volatility is extremely high, causing 1 minus the product to result in a negative value, the time confidence parameter is truncated to 0 rather than becoming negative, ensuring that the value range of the time confidence parameter always remains within the closed interval of 0 to 1. The time confidence parameter output in step S12 provides a time dimension constraint for the coupled product operation in step S13. Without step S12, the pixel-level dynamic confidence tensor generated in step S13 would be determined solely by the spatial confidence parameter. All pixels in the region near the center coordinates of the light source's brightness energy would receive high confidence, failing to distinguish between blood vessel pixels with stable grayscale gradients in the central region and pixels in the central region occasionally obscured by rapidly passing eyelash projections. Step S12 introduces temporal discriminative capabilities, complementing the spatial discriminative capabilities of step S11, enabling the subsequent pixel-level dynamic confidence tensor to simultaneously resist both spatial edge interference and temporal flicker interference through dual filtering.

[0063] Step S13: Perform coupled multiplication operation on the spatial confidence parameter and temporal confidence parameter of the target pixel to obtain the pixel-level dynamic confidence value of the target pixel, and assemble the pixel-level dynamic confidence values ​​of all target pixels in the current processing frame into the pixel-level dynamic confidence tensor of the current processing frame.

[0064] In step S13, the spatial confidence parameter of the target pixel obtained in step S11 and the temporal confidence parameter of the target pixel obtained in step S12 are subjected to a coupled product operation. The formula for the coupled product operation is as follows:

[0065] ;

[0066] In the formula, This represents the pixel-level dynamic confidence value of the target pixel. This represents the third scaling factor. The third scaling factor is a positive real number used to scale the coupled product result to match the numerical range of the confidence threshold set in step S21. The determination of the third scaling factor needs to consider the expected working range of the confidence threshold set in step S21. The third scaling factor should ensure that the pixel-level dynamic confidence value of the effective blood vessel region pixels falls within a range slightly higher than the set confidence threshold, while the pixel-level dynamic confidence value of the interference region pixels falls within a range lower than the set confidence threshold. Since the spatial confidence parameter and the temporal confidence parameter both have a closed range of 0 to 1, their product also falls within the range of 0 to 1. Once the confidence threshold set in step S21 is determined, the third scaling factor linearly scales the product result to maintain a suitable margin between the product value of the blood vessel region pixels and the set confidence threshold. In practice, the median of the product of the spatial confidence parameter and the temporal confidence parameter of the effective blood vessel region pixels in the calibration sequence is first used as the median of the effective product. Then, the set confidence threshold is divided by the median of the effective product and multiplied by a safety margin coefficient slightly greater than 1 to obtain the value of the third scaling factor. The safety margin coefficient ranges from 1.1 to 1.3, and its function is to ensure that the pixel-level dynamic confidence value of most pixels in the effective blood vessel region is higher than the set confidence threshold and will not be mistakenly excluded due to statistical fluctuations. The above coupled product operation is performed on all pixels in the current processing frame one by one. The pixel-level dynamic confidence values ​​of all pixels are arranged according to the row and column spatial position of each pixel in the current processing frame, generating a pixel-level dynamic confidence tensor of the current processing frame with the same row and column dimensions as the current processing frame. The position of each element in the pixel-level dynamic confidence tensor corresponds one-to-one with the position of the corresponding pixel in the current processing frame, and the value of each element represents the reliability of the corresponding pixel as a registration feature in the current frame and at the current time.

[0067] The mathematical structure of the coupled product operation determines that the pixel-level dynamic confidence value is only in the high-value range when both the spatial and temporal confidence parameters of the target pixel are high; when either the spatial or temporal confidence parameter is low, the pixel-level dynamic confidence value is suppressed to the low-value range. This product coupling characteristic produces a better distinguishing effect than comparing the two parameters separately with independent thresholds and then taking a logical AND, because the product operation preserves the continuous gradient information of the parameter values ​​instead of binarizing them. This allows pixels in "medium spatial location and medium temporal stability" to obtain medium confidence values, providing continuous weight gradients for the weighted least squares fitting in step S21 and avoiding matching jumps caused by hard threshold truncation at weight boundaries. Step S13 merges the spatial physical constraints of step S11 and the temporal statistical constraints of step S12 into a single tensor representation, so that the subsequent step S21 only needs to rely on a single tensor to make the decision of "retaining effective pixels and eliminating interfering pixels", without needing to embed a separate eyelash detection module or blood vessel segmentation module in the registration process.

[0068] Step S14: Perform block matching on the current processing frame and the previous frame to extract the mean and standard deviation of the displacement magnitude. Superimpose the mean and standard deviation of the displacement magnitude to output the global motion divergence parameter of the current processing frame. Compare the global motion divergence parameter with the first state boundary threshold and the second state boundary threshold to classify the current processing frame into the first state interval, the second state interval, or the third state interval as the motion state interval determination result.

[0069] See Figure 3 This is a schematic diagram of the block matching local displacement vector motion state interval division provided in an embodiment of this application, wherein... Figure 3 (a) illustrates the distribution of local displacement vectors generated after the current processing frame is divided into multiple non-overlapping rectangular sub-blocks. Figure 3 (b) illustrates three motion state intervals obtained by dividing the global motion divergence parameter using a dual threshold, presenting the range of global motion divergence parameters corresponding to different eye movement states. In step S14, block matching is performed on the current processing frame and the previous frame to obtain inter-frame motion information. The block matching operation divides the current processing frame into multiple non-overlapping rectangular sub-blocks. For each rectangular sub-block, the matching position with the smallest grayscale difference in the previous frame is searched. The coordinate difference between the matching position of each rectangular sub-block and its original position in the current processing frame constitutes a local displacement vector. Figure 3(a) is an example of the local displacement vector distribution obtained after dividing the current processing frame into rectangular sub-blocks. The magnitude of each local displacement vector in all rectangular sub-blocks is calculated. The arithmetic mean of the magnitudes of all local displacement vectors is calculated to obtain the mean displacement magnitude, and the standard deviation of the magnitudes of all local displacement vectors is calculated to obtain the standard deviation of the displacement magnitude. The mean displacement magnitude and the standard deviation of the displacement magnitude are added together to obtain the global motion divergence parameter of the current processing frame. The physical meaning of the global motion divergence parameter is a comprehensive representation of the overall motion amplitude and motion unevenness of the current processing frame relative to the previous frame: the mean displacement magnitude reflects the magnitude of the overall translation between frames, and the standard deviation of the displacement magnitude reflects the degree of inconsistency in the motion direction and amplitude of different regions between frames. When the eyeball is in a sub-pixel micro-tremor state, the magnitudes of the local displacement vectors of all sub-blocks are extremely small and similar to each other. The mean displacement magnitude and the standard deviation of the displacement magnitude are both extremely small, and the global motion divergence parameter is in a low value range, corresponding to... Figure 3 (b) The micro-tremor interval on the global motion divergence parameter axis is less than the first state boundary threshold; when the eyeball undergoes large-amplitude saccades, all sub-blocks exhibit large-amplitude displacements in the same direction, with a large mean displacement magnitude and a relatively controllable standard deviation of displacement magnitude due to the consistent direction. The global motion divergence parameter is in the medium-to-high range, corresponding to... Figure 3 (b) The global motion divergence parameter lies on the eye-tracking interval between the first and second state boundary thresholds. When the eyelids are fully closed, the grayscale matrix of the current frame undergoes a drastic global grayscale change due to occlusion compared to the previous frame. The local displacement vectors generated by the block matching operation exhibit random dispersion characteristics in both direction and amplitude. The mean and standard deviation of the displacement magnitudes are extremely high, and the global motion divergence parameter is in an extremely high range. Figure 3 (b) The blink interval on the global motion divergence parameter axis that is greater than the second state boundary threshold.

[0070] Based on the global motion divergence parameters of the current processing frame, the global motion divergence parameters are compared with the first state boundary threshold and the second state boundary threshold: when the value of the global motion divergence parameter is less than the first state boundary threshold, the current processing frame is determined to be in the first state interval, which represents the subpixel micro-tremor state; when the value of the global motion divergence parameter is greater than or equal to the first state boundary threshold and less than the second state boundary threshold, the current processing frame is determined to be in the second state interval, which represents the large-amplitude eye movement transition state; when the value of the global motion divergence parameter is greater than or equal to the second state boundary threshold, the current processing frame is determined to be in the third state interval, which represents the blink lock state. Since the inter-frame global motion amplitude represented by the first state interval, the second state interval, and the third state interval exhibits a progressively increasing physical characteristic, the first state boundary threshold is strictly less than the second state boundary threshold during the numerical setting and dynamic adjustment process. Both the first state boundary threshold and the second state boundary threshold are dynamic variables, and the initial values ​​of the first state boundary threshold and the second state boundary threshold, as well as the subsequent dynamic adjustment method, are determined by step S24. The initial value of the first state boundary threshold needs to consider two factors: the inherent jitter noise level of the imaging system and the typical amplitude of physiological eye tremors. If the initial value of the first state boundary threshold is too small and lower than the block matching error level caused by the inherent jitter noise of the imaging system, a large number of frames in the normal tremor state will be incorrectly classified into the second state interval and will not be able to participate in the fine topological registration mode of step S21, resulting in a sharp reduction in the number of frames available for high-quality fusion. If the initial value of the first state boundary threshold is too large, some frames that have already shown obvious eye movements will be incorrectly classified into the first state interval and enter the fine registration process of step S21, introducing registration errors. In actual determination, calibration sequences are collected under the stable gaze state of the test subject, and the distribution of global motion divergence parameters between adjacent frames in the calibration sequence is statistically analyzed. The values ​​in the range of the 95th to 99th quantiles at the upper end of the distribution are taken as the initial value of the first state boundary threshold. The initial value of the second-state boundary threshold needs to consider the prerequisite that block matching operations in large-amplitude eye-tracking frames can still generate meaningful local displacement vectors: if the initial value of the second-state boundary threshold is too small, some frames with obvious translation but still recognizable vascular textures will be classified into the third-state interval and directly discarded, wasting information that can be used for macroscopic position tracking; if the initial value of the second-state boundary threshold is too large, blinking frames will be classified into the second-state interval, and in step S22, an attempt will be made to extract the main vascular intersection anchor points, but there is no vascular information in the blinking frames, resulting in anchor point extraction failure or the extraction of false anchor points. In actual determination, the block matching divergence level caused by the sudden change in the grayscale matrix between adjacent frames when blinking occurs is used as the reference lower limit.When a blink occurs, the eyelid closes from top to bottom, causing a large area in the current processing frame to be covered by skin texture, resulting in a global grayscale difference compared to the fundus tissue image of the previous frame. At this time, in the block matching operation, each rectangular sub-block cannot find a region with a matching grayscale structure in the previous frame. The direction and amplitude of the local displacement vector exhibit an irregular random distribution, and the global motion divergence parameter jumps to a high value range due to a sharp increase in both the mean and standard deviation of the displacement magnitude. In a typical laser speckle fundus imaging system, the global motion divergence parameter of a complete blink frame usually exceeds 5 pixels and even reaches more than 10 pixels. The second state boundary threshold should be set between the upper limit of the divergence of large-amplitude eye-tracking frames and the lower limit of the divergence of blink frames. The global motion divergence parameter of large-amplitude eye-tracking frames is usually in the range of 1 to 4 pixels. Therefore, the initial value of the second state boundary threshold should be in the range of 2 to 5 pixels, so that large-amplitude eye-tracking frames are classified into the second state interval for macroscopic position tracking rather than being discarded as blink frames. For example, the initial value of the second state boundary threshold is 3 pixels.

[0071] The motion state interval determination result output in step S14 directly determines which of the three processing paths (steps S21, S22, and S23) in step S20 will be triggered. Without step S14, all frames in step S20 would indiscriminately enter the same registration path, and the incorrect matching results of large eye-tracking frames and blinking frames would contaminate the dynamically synthesized reference frame through the fusion update operation in step S24, leading to a continuous degradation in the image quality of the dynamically synthesized reference frame. Step S14 transforms the continuous values ​​of inter-frame motion into discrete state labels through a dual-threshold mechanism, providing clear branch triggering conditions for step S20. This ensures that the three functions—fine registration in step S21, macroscopic positioning in step S22, and trajectory deduction in step S23—are logically isolated from each other and do not overlap. Simultaneously, the first and second state boundary thresholds are set as dynamic variables rather than fixed constants, allowing feedback information on the reference frame quality in step S24 to be transmitted back to step S14, achieving a closed-loop adjustment capability where the motion state determination conditions dynamically adjust according to the actual acquisition quality.

[0072] Step S10, through the product coupling structure of spatial and temporal confidence parameters, unifies two fundamentally different image degradation mechanisms in laser speckle fundus images—eyelash occlusion interference and ocular bowl-shaped distortion—into a pure numerical computation framework that does not require image semantic understanding. This allows each pixel to obtain a quantitative scale representing the registration reliability without relying on any vascular morphology recognition or edge detection preprocessing. The spatial confidence parameter encodes the radially increasing physical law of distortion caused by the ocular bowl-shaped structure into a numerical mapping of exponential decay, while the temporal confidence parameter encodes the dynamic flicker characteristics of eyelash projection into a statistical deviation of gray-level gradient volatility. The pixel-level dynamic confidence tensor formed by the product operation of these two parameters numerically responds to both spatial location degradation and temporal stability degradation. The global motion divergence parameter takes into account the spatial non-uniformity of inter-frame motion rather than just calculating the average displacement of the entire frame, enabling the local differential motion caused by eye rotation and the global random dispersion motion caused by blinking to be distinguished into different state intervals. The combined output of the pixel-level dynamic confidence tensor and the motion state interval determination result enables step S20 to implement differentiated processing strategies at the frame level while simultaneously implementing differentiated weight allocation at the pixel level. This dual-layer architecture of frame-level coarse classification and pixel-level fine weighting breaks the traditional registration method's approach of "performing the same registration process uniformly for each frame." Under the premise of only increasing computational resource consumption by one gray-scale centroid calculation and one time sliding window statistical operation, it blocks the propagation path of interfering pixels and inferior frames to subsequent fusion stages from the source.

[0073] Step S20: Initialize the dynamic synthesis reference frame with the first frame of the laser speckle fundus image sequence. Trigger the corresponding registration processing path according to the motion state interval to which the current processing frame belongs. When the current processing frame is in the first state interval, perform fine topological registration and output the registered current frame. Calculate the fusion ratio based on the pixel-level dynamic confidence tensor. Perform weighted fusion update of the registered current frame and the dynamic synthesis reference frame according to the fusion ratio to obtain the updated dynamic synthesis reference frame. Adjust the first state boundary threshold and the second state boundary threshold based on the updated dynamic synthesis reference frame.

[0074] After step S10 outputs two types of information: pixel-by-pixel reliability measurement and frame-by-frame motion state label, these two types of information need to be transformed into specific differentiated processing actions at the registration execution level. This allows frames in the sub-pixel micro-tremor state to participate in deformation alignment and fusion with pixel-level precision, frames in the large eye-tracking state to only provide macro-position tracking without participating in pixel-level fusion, and frames in the blink-locked state to be completely isolated. At the same time, a dynamically synthesized reference frame that continuously evolves with the sequence processing process needs to be constructed to replace the manually selected fixed reference frame. The direction of change in the quality of the dynamically synthesized reference frame can inversely constrain the value of the motion state boundary threshold in step S14, forming a complete closed-loop adjustment structure from state determination to registration execution to quality evaluation.

[0075] Further, see Figure 4 Step S20 includes:

[0076] Step S21: When the motion state interval determination result indicates that the current processing frame is in the first state interval, the fine topological registration mode is triggered. Pixels with pixel-level dynamic confidence values ​​greater than the set confidence threshold are extracted from the pixel-level dynamic confidence tensor of the current processing frame to form a high-confidence pixel set. The dynamic synthesis reference frame is initialized with the first frame of the laser speckle fundus image sequence. The non-rigid deformation field between the current processing frame and the dynamic synthesis reference frame is fitted using the high-confidence pixel set. The current processing frame is resampled according to the non-rigid deformation field to output the registered current frame.

[0077] In step S21, when the motion state interval determination result output in step S14 indicates that the current processing frame is in the first state interval, the fine topology registration mode is triggered. From the pixel-level dynamic confidence tensor of the current processing frame obtained in step S13, the pixel-level dynamic confidence value is checked pixel by pixel to see if it is greater than the set confidence threshold. All pixels with pixel-level dynamic confidence values ​​greater than the set confidence threshold are collected into the high-confidence pixel set of the current processing frame. The determination of the confidence threshold needs to consider the degree of separation between the pixel-level dynamic confidence value distribution of pixels in the effective blood vessel region and the pixel-level dynamic confidence value distribution of pixels in the interference region: if the confidence threshold is set too low, some pixels in the interference region whose pixel-level dynamic confidence values ​​are not sufficiently suppressed will be included in the high-confidence pixel set. The false gray-level gradient information carried by these pixels will participate in the subsequent deformation field fitting and introduce erroneous deformation components; if the confidence threshold is set too high, only a very small number of pixels located at the center of the light source brightness energy center coordinates and with extremely low gray-level gradient fluctuation rate will be retained. The number of pixels in the high-confidence pixel set is insufficient to constrain the fitting degree of freedom of the non-rigid deformation field, resulting in overfitting or oscillation of the deformation field in the sparse pixel region. In actual determination, the cumulative distribution of all pixel-level dynamic confidence values ​​in the pixel-level dynamic confidence tensor output by statistical step S13 is used. The quantile value that makes the number of pixels at the high end of the cumulative distribution account for 30% to 50% of the total number of pixels in the current processing frame is selected as the set confidence threshold, so that the set of high-confidence pixels not only excludes obvious interference pixels, but also retains enough vascular region pixels for fitting.

[0078] The dynamic synthesis reference frame is initialized with the first frame of the laser speckle fundus image sequence when the laser speckle fundus image sequence processing starts. This means the grayscale matrix of the first frame is directly assigned as the initial grayscale matrix of the dynamic synthesis reference frame, and the number of fused frames is initialized to 1. The number of fused frames represents the cumulative number of frames that have participated in the weighted fusion of the dynamic synthesis reference frame up to the current time, and is used to adjust the denominator of the fusion weights in step S24. The first frame is used as the initialization object instead of intermediate frames in the sequence because the sequence processing proceeds sequentially from front to back, and the first frame can be obtained at the start of processing without waiting for subsequent frames to arrive. Based on the high-confidence pixel set, the gray-level gradient direction and magnitude in the current processing frame are extracted for each pixel in the high-confidence pixel set. In the dynamically synthesized reference frame, matching pixels with similar gray-level gradient direction and magnitude are searched within a search neighborhood centered on the corresponding pixel. The pixel with the smallest weighted sum of the differences in gray-level gradient direction and magnitude within the search neighborhood is determined as the matching pixel. This yields the displacement pair between the coordinates of each pixel in the high-confidence pixel set in the current processing frame and the coordinates of its corresponding matching pixel in the dynamically synthesized reference frame. The size of the search neighborhood is determined by the average displacement magnitude in step S14, with twice the average displacement magnitude used as the upper limit of the search neighborhood radius, ensuring that the search range covers the actual inter-frame displacement without introducing excessive search overhead.

[0079] Based on all displacement pairs, a weighted least squares fitting method is used to calculate the non-rigid deformation field between the current processing frame and the dynamically synthesized reference frame. The non-rigid deformation field characterizes the two-dimensional spatial offset required for each pixel in the current processing frame to align with the dynamically synthesized reference frame. The non-rigid deformation field is parameterized using a thin-plate spline interpolation model, a method in the field of scattered data interpolation used to construct smooth two-dimensional mappings. Known control point pairs are used as constraints to generate continuous deformation surfaces between the control points that satisfy the minimum bending energy criterion. During the weighted least squares fitting process, the fitting weight of each displacement pair is equal to the pixel-level dynamic confidence value of the corresponding pixel in the pixel-level dynamic confidence tensor. This results in displacement pairs with high confidence having a greater constraint on the non-rigid deformation field, and displacement pairs with low confidence having a smaller constraint. Instead of assigning equal weights to all displacement pairs, pixel-level dynamic confidence values ​​are used as fitting weights. This ensures that even if a small number of medium-confidence pixels near the boundaries are mixed into the high-confidence pixel set, the potential errors carried by these pixels are suppressed by low weights during the fitting process and do not dominate the direction of the deformation field. Based on the non-rigid deformation field, bilinear interpolation resampling is performed on each pixel in the current processing frame according to the two-dimensional offset given by the non-rigid deformation field to obtain the registered current frame aligned with the dynamically synthesized reference frame at the pixel level. A weighted mean displacement vector is calculated based on the magnitude and direction of all displacement pairs in the high-confidence pixel set. Using the pixel-level dynamic confidence value as weight, the whole-frame weighted rigid offset vector of the current processing frame relative to the dynamically synthesized reference frame is obtained. To ensure that the motion trajectory history sequence stores a uniform data type of rigid offset vector across all motion state intervals, a unified rigid offset vector for the entire frame is extracted from the weighted rigid offset vector of the entire frame. This is achieved by weighting all displacement pairs in the high-confidence pixel set according to their pixel-level dynamic confidence values, calculating the weighted median of all displacement pairs in the row direction and the weighted median in the column direction. The resulting two-dimensional vector, composed of the row and column weighted medians, is the unified rigid offset vector for the current processing frame. This unified rigid offset vector is recorded in the motion trajectory history sequence at the time position corresponding to the current processing frame, ensuring that the motion trajectory history sequence is updated frame-by-frame across multiple consecutive first state interval frames. This provides a temporally adjacent offset vector reference for step S23 in subsequent possible third state interval frames. Without step S21, high-quality frames in the first state interval would lose the channel to enter pixel-level deformation alignment, the dynamically synthesized reference frame would not receive new high-quality pixel information and would remain at the first frame level. The subsequent fusion update mechanism in step S24 would lose its input source, and the entire closed-loop structure would fail to start. Step S21 transforms the pixel-level dynamic confidence tensor output in step S13 from a passive metric into a weighting factor that actively participates in the fitting calculation, so that the reliability assessment results in step S10 can be substantially used in the registration execution stage of step S20.

[0080] Step S22: When the motion state interval determination result indicates that the current processing frame is in the second state interval, the rigid anchor point suspension mode is triggered. The region where the pixel-level dynamic confidence value is at a local maximum is located in the pixel-level dynamic confidence tensor of the current processing frame. Feature points where the gray-level gradient direction converges are extracted from the local maximum region as the main blood vessel intersection anchor points. The median of the coordinate offset between the main blood vessel intersection anchor points and the dynamic synthesis reference frame is calculated. The whole frame coarse rigid offset vector of the current processing frame is output. The whole frame coarse rigid offset vector is recorded as the whole frame unified rigid offset vector of the current processing frame in the motion trajectory history sequence and the current processing frame is prohibited from participating in pixel-level fusion.

[0081] In step S22, when the motion state interval determination result output in step S14 indicates that the current processing frame is in the second state interval, the rigid anchor point suspension mode is triggered. (See also...) Figure 5 This is a schematic diagram of the extraction of the main blood vessel intersection anchor points provided in the embodiments of this application. Figure 5 The image illustrates the structure of the optic disc and the distribution of the main blood vessels originating from the optic disc in a fundus image. Figure 5 The circular black dots represent the main vessel intersections or bifurcations. In the pixel-level dynamic confidence tensor of the current processing frame, local maxima detection is performed: a preset neighborhood window is taken centered on each pixel. When the pixel-level dynamic confidence value of the center pixel is greater than the pixel-level dynamic confidence values ​​of all adjacent pixels within the preset neighborhood window, the location of the center pixel is determined to be a local maximum region of the pixel-level dynamic confidence value. The side length of the preset neighborhood window needs to consider the typical width and bifurcation distance of the main vessels in the fundus: if the side length of the preset neighborhood window is too small, there will be too many and too dense local maximum regions, and many single vessel segments not at vessel intersections will also be misclassified as local maximum regions; if the side length of the preset neighborhood window is too large, only a very few maximum regions will be retained, which may miss important vessel intersection structures located at the edge of the optic disc. In practice, the side length of the preset neighborhood window is set to 3 to 5 times the typical width of the main retinal vessels, ensuring that each local maximum region corresponds to a spatially independent and prominent vascular structural feature with a high confidence level. At the resolution of common laser speckle retinal imaging systems, the width of the main retinal vessels (i.e., the primary arteries and veins originating from the optic disc) in the image is typically 4 to 7 pixels. Multiplying this width by 3 to 5 yields a range of 12 to 35 pixels, and the middle segment of this range is taken as the side length of the preset neighborhood window. For example, the side length of the preset neighborhood window is set to 15 to 25 pixels.

[0082] Within the local maximum region, the gray-level gradient direction is calculated for each pixel. Feature points whose gray-level gradient directions spatially converge towards the center are extracted as anchor points for the main blood vessel intersections.Figure 5 The area marked by the dashed circle represents a typical anchor point where major blood vessels intersect. The arrows indicate that the grayscale gradient direction of each pixel within this area all points towards the intersection center. The physical meaning of this convergence of grayscale gradient directions is that the grayscale at the intersection or bifurcation of blood vessels decreases from the vessel wall towards the center, and the grayscale gradient directions of multiple blood vessels all point towards the intersection center. Figure 5The optic disc, as the largest vascular convergence center, also exhibits the characteristic of all main vascular grayscale gradient directions converging towards the center of the optic disc. This feature remains recognizable even when the image is blurred due to significant eye movement, because significant motion blur only weakens the magnitude of the grayscale gradient without changing the convergence topology. Based on the coordinates of each main vascular convergence anchor point extracted in the current processing frame and the corresponding matching coordinates in the dynamic synthesis reference frame, within the dynamic synthesis reference frame, centered on the coordinates of the anchor point in the current processing frame and with a search radius of three times the average displacement magnitude in step S14, the corresponding vascular convergence feature points in the dynamic synthesis reference frame are located using the same grayscale gradient direction convergence detection method as in the current processing frame. The coordinates of the feature point with the highest degree of grayscale gradient convergence are determined as the matching coordinates. The coordinate offset of each main vascular convergence anchor point is calculated, and the median of the coordinate offsets of all main vascular convergence anchor points is taken to obtain the overall coarse rigid offset vector of the current processing frame. The median calculation is used instead of the mean calculation because some anchor points in large-scale eye-tracking frames may have outlier offset values ​​due to local occlusion or blurring. The median calculation is naturally robust to outliers and is not dominated by a single erroneous anchor point. The whole-frame coarse rigid offset vector is used as the whole-frame unified rigid offset vector of the current processing frame, and the whole-frame unified rigid offset vector is recorded in the motion trajectory history sequence at the time position corresponding to the current processing frame. The motion trajectory history sequence is a set of frame-level displacement vectors that store the whole-frame unified rigid offset vectors of each frame in chronological order. The whole-frame unified rigid offset vector of the first state interval frame in the motion trajectory history sequence is extracted from the displacement pairs of the non-rigid deformation field in step S21, and the whole-frame unified rigid offset vector of the second state interval frame is calculated from the median of the coordinate offsets of the main vessel intersection anchor points in step S22. Both types of whole-frame unified rigid offset vectors are stored in the motion trajectory history sequence using the same data format and are used for trajectory extrapolation calculation in step S23. The current processing frame is prohibited from participating in pixel-level fusion, that is, it does not enter the weighted fusion update process in step S24. Step S22 isolates large-amplitude eye-tracking frames from the pixel-level registration process, allowing only identifiable coarse topological information from the intersection of the optic disc and main blood vessels to participate in macroscopic translation estimation. This prevents false detail features in the blurred image from being mixed into the dynamically synthesized reference frame through the fusion operation in step S24. Without step S22, large-amplitude eye-tracking frames would be discarded as blink frames in step S23, causing positional breaks in the motion trajectory history sequence during large-amplitude eye movements. The trajectory extrapolation in step S23 after the blink would then suffer from extrapolation bias due to the lack of true positional information during large-amplitude eye movements. Step S22 establishes an intermediate processing layer between the fine registration in step S21 and the complete discarding in step S23, ensuring the temporal continuity of the motion trajectory history sequence across all motion states.

[0083] Step S23: When the motion state interval determination result indicates that the current processing frame is in the third state interval, the trajectory extrapolation skip mode is triggered, the image information of the current processing frame is discarded, and the predicted position coordinates of the current processing frame are estimated by linear extrapolation based on the unified rigid offset vector of the whole frame recorded in the motion trajectory history sequence. The predicted position coordinates are then added to the motion trajectory history sequence to maintain the temporal continuity of the motion trajectory history sequence.

[0084] In step S23, when the motion state interval determination result output in step S14 indicates that the current processing frame is in the third state interval, the trajectory extrapolation skip mode is triggered. The grayscale matrix information of the current processing frame is discarded, and no feature extraction, matching, or deformation calculation operations are performed on the current processing frame. Based on the whole-frame unified rigid offset vector recorded in the motion trajectory history sequence, the two consecutive whole-frame unified rigid offset vectors that are closest to the current processing frame in time in the motion trajectory history sequence are taken, and the difference between the two consecutive whole-frame unified rigid offset vectors is calculated as the inter-frame velocity estimate. The predicted position coordinates of the current processing frame are obtained by adding the whole-frame unified rigid offset vector closest to the current processing frame to the inter-frame velocity estimate. Since the first state interval frame and the second state interval frame in the motion trajectory history sequence both store the whole-frame unified rigid offset vector, regardless of how the motion state interval determination results before the current processing frame are distributed, there are continuous offset vector records in the motion trajectory history sequence that can be used for linear extrapolation. This avoids the problem of lacking extrapolation data when encountering the third state interval after being in the first state interval for a long time, which is only based on the stored data of a single motion state interval. The aforementioned linear extrapolation method is based on the physical assumption that eye movement speed remains stable over a short period. Since a blink typically lasts only a few frames, linear extrapolation is sufficient to provide reasonable position prediction accuracy. The predicted position coordinates are then added to the time position of the corresponding current processing frame in the motion trajectory history sequence, ensuring the temporal index of the motion trajectory history sequence remains continuous and without gaps. The complete isolation of the blink frame in step S23 prevents global grayscale abrupt changes caused by eyelid occlusion from being misinterpreted as large spatial displacements. If step S23 is omitted and the blink frame is sent to the anchor point extraction process in step S22, the blink frame contains no vascular structure information. The anchor point extraction algorithm will locate false convergence points in the random grayscale texture and generate incorrect offset vectors, contaminating the motion trajectory history sequence. Step S23 maintains the continuity of the motion trajectory history sequence during the blink interval through trajectory extrapolation, enabling the first first-state interval frame or the second-state interval frame immediately following the blink to quickly recover the spatial correspondence with the dynamically synthesized reference frame based on continuous trajectory information.

[0085] Step S24: Calculate the fusion ratio based on the dynamic confidence value in the pixel-level dynamic confidence tensor. Perform weighted fusion update on the registered current frame and the dynamic synthesis reference frame according to the fusion ratio to synthesize the updated dynamic synthesis reference frame. Calculate the difference between the image information entropy value of the updated dynamic synthesis reference frame and the image information entropy value of the dynamic synthesis reference frame before the update to obtain the entropy change. Adjust the first state boundary threshold and the second state boundary threshold according to the positive and negative direction and absolute value of the entropy change.

[0086] In step S24, when the registered current frame is output in step S21, a weighted fusion update of the dynamically synthesized reference frame is performed. For each pixel in the registered current frame, the second fusion ratio and the first fusion ratio of the corresponding pixel are calculated. The calculation of the second fusion ratio is divided into two steps: First, the candidate value of the second fusion ratio is calculated. The candidate value of the second fusion ratio is equal to the pixel-level dynamic confidence value multiplied by the reciprocal of the number of fused frames and then multiplied by the fourth scaling factor. Wherein, the pixel-level dynamic confidence value is the value of the corresponding pixel position in the pixel-level dynamic confidence tensor of the current processing frame output in step S13, the number of fused frames is the cumulative number of times the weighted fusion update of step S24 has been performed up to the current time plus 1, and the fourth scaling factor is a positive real number greater than 0 and less than or equal to 1. Then, an upper limit truncation operation is performed on the candidate value of the second fusion ratio, and the smaller value between the candidate value of the second fusion ratio and 1 is taken as the final value of the second fusion ratio. The upper limit truncation operation ensures that the value of the second fusion ratio is always within the closed interval of 0 to 1. The determination of the fourth scaling factor needs to consider the balance between the speed at which the fusion process incorporates new frame information and the degree to which it preserves existing reference frame information. If the fourth scaling factor is close to the upper limit of its range, the contribution of new frame pixels to the dynamically synthesized reference frame is relatively high in each fusion, and the frame information fused in the early stages is rapidly diluted. The speed at which the signal-to-noise ratio (SNR) of the dynamically synthesized reference frame improves is limited by the SNR of a single frame rather than by the accumulation of multiple frames. If the fourth scaling factor is close to the lower limit of its range, the contribution of new frame pixels to the dynamically synthesized reference frame is extremely low, and the image content of the dynamically synthesized reference frame remains near the first frame for a long time, failing to absorb new regional information brought about by micro-vibration displacement in subsequent frames. In practice, the fourth scaling factor is selected with the goal of maintaining a monotonically increasing SNR gain in the dynamically synthesized reference frame before it saturates after accumulating 10 to 20 frames of fusion. The fourth scaling factor controls the upper limit of the weight ratio of new frame information in the dynamically synthesized reference frame in each fusion. When the fourth scaling factor equals 1, in the early stages when the number of merged frames is small, the weights of the new frame and the existing reference frame are roughly equal, and the fusion effect is close to an equal-weighted average. As the number of merged frames increases, the weight of the new frame gradually decreases due to being divided by the number of merged frames, and the dynamically synthesized reference frame tends to stabilize. When the fourth scaling factor is less than 1, the injection speed of new frame information is slowed down, making the dynamically synthesized reference frame more inclined to maintain the stability of existing information. Since the value range of the fourth scaling factor is constrained to be greater than zero and less than or equal to 1, the candidate value of the second fusion ratio mathematically does not exceed the pixel-level dynamic confidence value multiplied by the reciprocal of the number of merged frames. The pixel-level dynamic confidence value ranges from 0 to 1 in a closed interval, and the number of merged frames is greater than or equal to 1. Therefore, the candidate value of the second fusion ratio does not exceed 1 under the condition that the number of merged frames is greater than or equal to 1. The upper limit truncation operation does not change the candidate value of the second fusion ratio under normal working conditions, but only serves as a numerical safety guarantee mechanism to prevent the second fusion ratio from exceeding the effective range under extreme parameter combinations.For laser speckle fundus imaging systems with frame rates of 60 to 80 frames per second, the ratio of vascular signal to speckle noise typically begins to show a significant slowdown in improvement after accumulating 8 to 15 frames of fusion. The fourth scaling factor should be selected such that each fusion within this frame range still results in an observable improvement in signal-to-noise ratio. Based on the above considerations, the value of the fourth scaling factor ranges from 0.5 to 1.0. When the speckle noise level of the imaging system is high, the value is closer to the upper limit to accelerate noise smoothing; when the speckle noise level is low, the value is closer to the lower limit to maintain the stability of the existing reference frame. For example, the fourth scaling factor is set to 0.8. The first fusion ratio equals 1 minus the second fusion ratio. Since the value of the second fusion ratio, after upper limit truncation, never exceeds 1, and the value of the first fusion ratio is always greater than or equal to zero, the weighting coefficients of the grayscale value of the current dynamically synthesized reference frame and the grayscale value of the registered current frame in the weighted fusion formula are both non-negative, ensuring the physical rationality of the weighted fusion operation.

[0087] In the updated dynamic composite reference frame, the grayscale value of each pixel is equal to the grayscale value of the corresponding pixel in the current dynamic composite reference frame multiplied by the first fusion ratio, plus the grayscale value of the corresponding pixel in the registered current frame multiplied by the second fusion ratio. Since the second fusion ratio includes a pixel-level dynamic confidence value as a factor, for pixels with high confidence, the new frame information contributes significantly to the dynamic composite reference frame; for pixels with low confidence, the new frame information contributes less, and the dynamic composite reference frame maintains its original grayscale value in low-confidence regions, unaffected by the new frame. The number of fused frames increments by 1 after each execution of step S24.

[0088] After completing the weighted fusion update of the dynamically synthesized reference frame, the image information entropy value of the updated dynamically synthesized reference frame is calculated. The image information entropy value is calculated using the Shannon information entropy formula on the gray-level histogram. The probability of occurrence of each gray level in the gray-level histogram of the updated dynamically synthesized reference frame is multiplied by the negative logarithm (base 2) of the corresponding probability, and then summed to obtain the updated reference frame entropy value. Simultaneously, the image information entropy value of the dynamic synthesized reference frame before the update is retained as the original reference frame entropy value. The image information entropy value characterizes the richness of the image's gray-level distribution. A higher image information entropy value indicates clear blood vessel textures and distinct gray-level levels, while a lower value indicates a blurred image or uniform noise contamination. The updated reference frame entropy value is compared with the original reference frame entropy value, and the change in entropy value is calculated by subtracting the original reference frame entropy value from the updated reference frame entropy value. When the entropy change is greater than 0, it indicates that the fusion of the current frame after registration has improved the image quality of the dynamically synthesized reference frame. In this case, the first state boundary threshold in step S14 is increased by the first adjustment step size, and the second state boundary threshold is increased by the second adjustment step size. When the entropy change is less than or equal to zero, it indicates that the fusion of the current frame after registration has failed to improve the image quality of the dynamically synthesized reference frame. In this case, the first state boundary threshold in step S14 is decreased by the first adjustment step size, and the second state boundary threshold is decreased by the second adjustment step size. The first adjustment step size is equal to the absolute value of the entropy change multiplied by the fifth scaling factor, and the second adjustment step size is equal to the first adjustment step size multiplied by the sixth scaling factor.

[0089] The determination of the fifth scaling factor needs to consider the trade-off between the sensitivity and stability of the threshold adjustment response: if the fifth scaling factor is too large, even a small fluctuation in entropy caused by a single fusion can lead to a significant jump in the first state boundary threshold, causing the motion state determination results of adjacent frames to frequently switch between the first and second state intervals, disrupting the temporal continuity of the registration process; if the fifth scaling factor is too small, the threshold response to changes in image quality is too slow, requiring multiple frames to accumulate before tightening the judgment conditions in the event of a sudden decrease in shooting coordination, during which a large number of edge-quality frames have already been mixed into the dynamically synthesized reference frame. In practice, the threshold adjustment amount resulting from the accumulation of five consecutive entropy changes in the same direction is selected with the goal of being within 10% to 20% of the initial value of the first state boundary threshold. The determination of the sixth scaling factor needs to consider the proportional relationship between the first state boundary threshold and the second state boundary threshold: if the sixth scaling factor is too small, the adjustment range of the second state boundary threshold is much smaller than that of the first state boundary threshold, causing the distance between the two thresholds to gradually shrink until they overlap after a long period of processing, resulting in the disappearance of the second state interval and the loss of the triggering opportunity for step S22; if the sixth scaling factor is too large, the adjustment range of the second state boundary threshold is much larger than that of the first state boundary threshold, causing the coverage of the second state interval to expand or contract excessively, thus losing the ability to accurately capture large eye-tracking frames. In actual determination, the ratio of the initial value of the second state boundary threshold to the initial value of the first state boundary threshold is used as a reference to keep the sixth scaling factor and the ratio of the initial values ​​of the two thresholds on the same order of magnitude. The reason for selecting this order of magnitude is that the first state boundary threshold and the second state boundary threshold need to maintain the stability of their relative distance during the adjustment process, that is, to always maintain a fixed multiple relationship between them, to avoid the situation where the growth rate of the first state boundary threshold is much faster or slower than that of the second state boundary threshold, resulting in excessive contraction or expansion of the second state interval.

[0090] When the entropy change is greater than 0, the first and second state boundary thresholds are increased, causing more frames to be identified as belonging to the first state interval and enter the fine registration path in step S21. The logic behind this adjustment is that improved fusion quality indicates that frames entering fine registration under the current judgment conditions are indeed high-quality frames, and the conditions can be appropriately relaxed to allow more frames to participate in fusion to accelerate signal-to-noise ratio accumulation. When the entropy change is less than or equal to 0, the first and second state boundary thresholds are decreased, tightening the judgment conditions and reducing the chance of subsequent frames entering step S21. The logic behind this adjustment is that no improvement in fusion quality indicates that frames entering fine registration may contain frames with poor registration quality, and the admission conditions need to be tightened to protect the dynamically synthesized reference frame from further degradation. The closed-loop feedback structure constructed in step S24 makes the first and second state boundary thresholds no longer static parameters that remain unchanged throughout the processing, but variables that are dynamically adjusted according to the actual acquisition conditions. When facing test subjects with high cooperation during filming, most frames are in a state of slight tremor, and the fusion quality continues to improve. The first state boundary threshold is gradually widened, allowing more frames to participate in fine registration, and the signal-to-noise ratio of the dynamically synthesized reference frame accumulates rapidly. When facing test subjects with low cooperation during filming, frequent eye movements cause fluctuations or even a decrease in fusion quality. The first state boundary threshold is gradually tightened, allowing only frames with extremely small tremor amplitudes to participate in fine registration, and the quality of the dynamically synthesized reference frame is strictly protected. This adaptive adjustment capability allows the same registration method to adapt to the differences in eye movement behavior of different test subjects without manual parameter modification. Step S24 extends the pixel-level dynamic confidence tensor output in step S13 from the matching weights in step S21 to fusion weights, so that the same set of reliability metrics plays a dual constraint function in the registration and fusion stages. In the registration stage, it constrains the direction of deformation field fitting, and in the fusion stage, it constrains the mixing ratio of new and old information. If step S24 is missing, the registered current frame output in step S21 will lack a channel to integrate with the dynamically synthesized reference frame. The dynamically synthesized reference frame will remain in the initial state of the first frame and will not be able to continuously absorb new information as the sequence progresses. The first state boundary threshold and the second state boundary threshold in step S14 will lose the driving source for dynamic adjustment and remain at a static initial value. The entire processing flow will degenerate into an open-loop execution mode without feedback.

[0091] Step S20 uses three mutually exclusive processing paths to separate frames of varying quality in the fundus image sequence: Step S21 performs non-rigid deformation field fitting with pixel-level confidence as weight to achieve fine alignment of micro-tremor frames; Step S22 uses the median offset of the main vascular intersection anchor points to achieve macroscopic position recording of large eye movement frames without contaminating the fusion result; Step S23 uses linear extrapolation to achieve trajectory continuation of blink frames without introducing false image information. Step S24 embeds pixel-level dynamic confidence values ​​into the fusion ratio formula, allowing the update rate of the dynamically synthesized reference frame at each pixel position to be independently controlled by the reliability of the corresponding pixel: high-confidence regions rapidly absorb new frame information, continuously clarifying vascular texture; low-confidence regions, due to low pixel-level dynamic confidence values, have a second fusion ratio approaching zero, almost completely excluding new frame information, thus preventing interference and noise from being solidified. This pixel-by-pixel differentiated fusion rate produces results that traditional equal-weighted averaging methods cannot achieve. The clarity of blood vessel texture in the dynamically synthesized reference frame increases monotonically with the number of fused frames, while the grayscale of the interference area remains in the initial frame state and is not covered by eyelash occlusion that may occur in subsequent frames. The closed-loop structure formed between the entropy feedback mechanism in step S24 and the dual threshold determination in step S14 establishes a causal relationship between the strictness of motion state determination and the actual quality performance of the dynamically synthesized reference frame. This breaks the open-loop limitation of traditional registration methods, which lack an information loop between the determination threshold and the registration result. This enables the entire registration system to self-correct based on the actual quality distribution of each acquired sequence without external manual intervention.

[0092] Example 2:

[0093] This embodiment, based on Embodiment 1, provides a dynamic fundus image registration system based on spatiotemporal coupling, such as... Figure 6 As shown, it includes:

[0094] Confidence Tensor Generation Module: Performs source brightness and energy center coordinate localization and time sliding window grayscale gradient fluctuation rate calculation on the current processing frame in the laser speckle fundus image sequence. Couples and multiplies the spatial confidence parameter driven by the source brightness and energy center coordinate with the time confidence parameter driven by the grayscale gradient fluctuation rate to generate the pixel-level dynamic confidence tensor of the current processing frame.

[0095] Motion state determination module: Performs block matching calculation on the current processing frame and the previous frame to obtain global motion divergence parameters, and compares the global motion divergence parameters with the first state boundary threshold and the second state boundary threshold to determine the motion state interval to which the current processing frame belongs;

[0096] Reference frame registration trigger module: Initializes the dynamic synthesis reference frame with the first frame of the laser speckle fundus image sequence, and triggers the corresponding registration processing path according to the motion state interval to which the current processing frame belongs;

[0097] Fine registration and fusion module: When the current processing frame is in the first state interval, fine topological registration is performed and the registered current frame is output. The fusion ratio is calculated based on the pixel-level dynamic confidence tensor. The registered current frame and the dynamic synthesis reference frame are weighted and fused according to the fusion ratio to obtain the updated dynamic synthesis reference frame. The first state boundary threshold and the second state boundary threshold are adjusted based on the updated dynamic synthesis reference frame.

[0098] Furthermore, in the confidence tensor generation module, the method for locating the coordinates of the light source brightness energy center includes:

[0099] Read the gray values ​​of all pixels in the current processing frame of the laser speckle fundus image sequence, and perform weighted centroid calculation on the gray values ​​of all pixels in the current processing frame to obtain the coordinates of the light source brightness energy center of the current processing frame.

[0100] The method for calculating the time-sliding window gray-scale gradient volatility includes:

[0101] Using the current processing frame as the last frame, the preceding consecutive frames are taken to form a time sliding window. For each target pixel within the time sliding window, the local gray-level gradient magnitude sequence is extracted. The standard deviation and mean of the local gray-level gradient magnitude sequence are calculated. The gray-level gradient fluctuation rate of the target pixel is obtained based on the ratio of the standard deviation to the mean of the gray-level gradient magnitude sequence.

[0102] The method for generating the pixel-level dynamic confidence tensor of the current processing frame includes:

[0103] The spatial confidence parameters of the target pixel are constructed based on the coordinates of the light source brightness energy center, and the temporal confidence parameters of the target pixel are constructed based on the gray-scale gradient fluctuation rate.

[0104] The spatial confidence parameter and temporal confidence parameter of the target pixel are coupled and multiplied to obtain the pixel-level dynamic confidence value of the target pixel. The pixel-level dynamic confidence values ​​of all target pixels in the current processing frame are then combined to form the pixel-level dynamic confidence tensor of the current processing frame.

[0105] Furthermore, in the motion state determination module, the calculation method for the global motion divergence parameter includes:

[0106] The current processing frame is divided into multiple non-overlapping rectangular sub-blocks. For each rectangular sub-block, the matching position with the smallest grayscale difference in the previous frame is searched to obtain the local displacement vector. The arithmetic mean of the magnitudes of all local displacement vectors is calculated to obtain the displacement magnitude mean. The standard deviation of the magnitudes of all local displacement vectors is calculated to obtain the displacement magnitude standard deviation. The displacement magnitude mean and the displacement magnitude standard deviation are added together to obtain the global motion divergence parameter.

[0107] The methods and systems of this application may be implemented in many ways. For example, they may be implemented by software, hardware, firmware, or any combination of software, hardware, and firmware. The above-described order of steps for the method is for illustrative purposes only, and the steps of the method of this application are not limited to the order specifically described above, unless otherwise specifically stated.

[0108] In addition, the parts of the technical solutions provided in the embodiments of this application that are consistent with the implementation principles of the corresponding technical solutions in the prior art have not been described in detail, so as to avoid excessive elaboration.

[0109] The specific embodiments described above further illustrate the purpose, technical solution, and beneficial effects of the present invention. It should be understood that the above descriptions are merely specific embodiments of the present invention and are not intended to limit the invention. Any modifications, equivalent substitutions, or improvements made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

[0110] It should be noted that the collection, storage, and use of all personal information involved in the technical solution of this invention must be carried out only after obtaining the explicit authorization and separate consent of the information subject. The processing of personal information-related data strictly complies with the requirements of currently effective national laws, regulations, and relevant standards and specifications. The collected personal information is limited to use within the specific purpose necessary to achieve the technical solution of this invention, and no processing beyond that purpose is performed. Necessary technical and management measures are taken to ensure the security of personal information.

Claims

1. A dynamic fundus image registration method based on spatiotemporal coupling, characterized in that, The method includes: For the current processing frame in the laser speckle fundus image sequence, the coordinates of the light source brightness and energy center are located and the gray-level gradient fluctuation rate of the time sliding window is calculated. The spatial confidence parameter driven by the light source brightness and energy center coordinates and the temporal confidence parameter driven by the gray-level gradient fluctuation rate are coupled and multiplied to generate the pixel-level dynamic confidence tensor of the current processing frame. The global motion divergence parameter is calculated and compared with the first state boundary threshold and the second state boundary threshold to determine the motion state interval to which the current processing frame belongs. The first frame of the laser speckle fundus image sequence is used to initialize the dynamic synthesis reference frame. The corresponding registration processing path is triggered according to the motion state interval to which the current processing frame belongs. This includes performing fine topological registration and outputting the registered current frame when the current processing frame is in the first state interval. The fusion ratio is calculated based on the pixel-level dynamic confidence tensor. The registered current frame and the dynamic synthesis reference frame are weighted and fused according to the fusion ratio to obtain the updated dynamic synthesis reference frame. The first state boundary threshold and the second state boundary threshold are adjusted based on the updated dynamic synthesis reference frame.

2. The dynamic fundus image registration method based on spatiotemporal coupling according to claim 1, characterized in that, The method for locating the center coordinates of the light source's brightness energy includes: Read the gray values ​​of all pixels in the current processing frame of the laser speckle fundus image sequence, and perform weighted centroid calculation on the gray values ​​of all pixels in the current processing frame to obtain the coordinates of the light source brightness energy center of the current processing frame.

3. The dynamic fundus image registration method based on spatiotemporal coupling according to claim 1, characterized in that, The method for calculating the time-sliding window gray-scale gradient volatility includes: Using the current processing frame as the last frame, the preceding consecutive frames are taken to form a time sliding window. For each target pixel within the time sliding window, the local gray-level gradient magnitude sequence is extracted. The standard deviation and mean of the local gray-level gradient magnitude sequence are calculated. The gray-level gradient fluctuation rate of the target pixel is obtained based on the ratio of the standard deviation to the mean of the gray-level gradient magnitude sequence.

4. The dynamic fundus image registration method based on spatiotemporal coupling according to claim 3, characterized in that, The method for generating the pixel-level dynamic confidence tensor of the current processing frame includes: The spatial confidence parameters of the target pixel are constructed based on the coordinates of the light source brightness energy center, and the temporal confidence parameters of the target pixel are constructed based on the gray-scale gradient fluctuation rate. The spatial confidence parameter and temporal confidence parameter of the target pixel are coupled and multiplied to obtain the pixel-level dynamic confidence value of the target pixel. The pixel-level dynamic confidence values ​​of all target pixels in the current processing frame are then combined to form the pixel-level dynamic confidence tensor of the current processing frame.

5. The dynamic fundus image registration method based on spatiotemporal coupling according to claim 4, characterized in that, The method for constructing the spatial confidence parameters of the target pixel based on the coordinates of the light source brightness energy center includes: The center offset distance of the target pixel is obtained by calculating the Euclidean distance between the target pixel and the center coordinates of the light source brightness energy in the current processing frame. The spatial confidence parameter of the target pixel is constructed based on the center offset distance.

6. The dynamic fundus image registration method based on spatiotemporal coupling according to claim 1, characterized in that, The calculation method for the global motion divergence parameter includes: The current processing frame is divided into multiple non-overlapping rectangular sub-blocks. For each rectangular sub-block, the matching position with the smallest grayscale difference in the previous frame is searched to obtain the local displacement vector. The arithmetic mean of the magnitudes of all local displacement vectors is calculated to obtain the displacement magnitude mean. The standard deviation of the magnitudes of all local displacement vectors is calculated to obtain the displacement magnitude standard deviation. The displacement magnitude mean and the displacement magnitude standard deviation are added together to obtain the global motion divergence parameter.

7. The dynamic fundus image registration method based on spatiotemporal coupling according to claim 6, characterized in that, The method for determining the motion state interval to which the current processing frame belongs includes: When the value of the global motion divergence parameter is less than the first state boundary threshold, the current processing frame is determined to be in the first state interval, which represents the subpixel micro-tremor state; when the value of the global motion divergence parameter is greater than or equal to the first state boundary threshold and less than the second state boundary threshold, the current processing frame is determined to be in the second state interval, which represents the large-amplitude eye movement transition state; when the value of the global motion divergence parameter is greater than or equal to the second state boundary threshold, the current processing frame is determined to be in the third state interval, which represents the blink lock state; wherein, the first state boundary threshold is less than the second state boundary threshold.

8. The dynamic fundus image registration method based on spatiotemporal coupling according to claim 1, characterized in that, The method for performing fine-grained topology registration includes: Pixels with pixel-level dynamic confidence values ​​greater than a set confidence threshold are extracted from the pixel-level dynamic confidence tensor of the current processing frame to form a high-confidence pixel set; the high-confidence pixel set is used to fit the non-rigid deformation field between the current processing frame and the dynamically synthesized reference frame, and the current processing frame is resampled according to the non-rigid deformation field to output the registered current frame.

9. The dynamic fundus image registration method based on spatiotemporal coupling according to claim 1, characterized in that, The step of triggering the corresponding registration processing path based on the motion state interval to which the current processing frame belongs also includes: When the motion state interval determination result indicates that the current processing frame is in the second state interval, the rigid anchor suspension mode is triggered; when the motion state interval determination result indicates that the current processing frame is in the third state interval, the trajectory deduction skip mode is triggered. The rigid anchor point suspension mode includes extracting the main vessel intersection anchor points from the pixel-level dynamic confidence tensor of the current processing frame, calculating the median of the coordinate offset between the main vessel intersection anchor points and the dynamic synthesis reference frame, outputting the whole frame coarse rigid offset vector of the current processing frame, recording the whole frame coarse rigid offset vector as the whole frame unified rigid offset vector of the current processing frame to the motion trajectory history sequence, and prohibiting the current processing frame from participating in pixel-level fusion.

10. A dynamic fundus image registration system based on spatiotemporal coupling, used to implement the dynamic fundus image registration method based on spatiotemporal coupling as described in any one of claims 1-9, characterized in that, The system includes: Confidence Tensor Generation Module: Performs source brightness and energy center coordinate localization and time sliding window grayscale gradient fluctuation rate calculation on the current processing frame in the laser speckle fundus image sequence. Couples and multiplies the spatial confidence parameter driven by the source brightness and energy center coordinate with the time confidence parameter driven by the grayscale gradient fluctuation rate to generate the pixel-level dynamic confidence tensor of the current processing frame. Motion state determination module: Performs block matching calculation on the current processing frame and the previous frame to obtain global motion divergence parameters, and compares the global motion divergence parameters with the first state boundary threshold and the second state boundary threshold to determine the motion state interval to which the current processing frame belongs; Reference frame registration trigger module: Initializes the dynamic synthesis reference frame with the first frame of the laser speckle fundus image sequence, and triggers the corresponding registration processing path according to the motion state interval to which the current processing frame belongs; Fine registration and fusion module: When the current processing frame is in the first state interval, fine topological registration is performed and the registered current frame is output. The fusion ratio is calculated based on the pixel-level dynamic confidence tensor. The registered current frame and the dynamic synthesis reference frame are weighted and fused according to the fusion ratio to obtain the updated dynamic synthesis reference frame. The first state boundary threshold and the second state boundary threshold are adjusted based on the updated dynamic synthesis reference frame.

Images