Hartmann wavefront sensor time super-resolution reconstruction method fusing physical constraints and based on pixel-level three-spline interpolation

Abstract

The invention relates to a Hartmann wavefront sensor time super-resolution reconstruction method fusing physical constraints and based on pixel-level three-spline interpolation, and belongs to the technical field of light beam wavefront measurement. Non-physical abrupt change caused by noise, sampling abnormity or pixel clamping is eliminated, and continuity and reliability of input data are ensured. On the basis, a time three-spline interpolation model is constructed at the pixel level, and multi-level physical constraints are introduced in the interpolation process, so that the reconstructed high-time-resolution light spot sequence is continuous and smooth in mathematics, and meanwhile, the reconstructed high-time-resolution light spot sequence more conforms to the actual wavefront evolution characteristic physically. According to the wavefront time super-resolution reconstruction method, the spatial-temporal correlation of the original light spot image can be fully utilized, and wavefront time super-resolution reconstruction with higher precision and higher robustness is realized.

Description

TECHNICAL FIELD

[0001] The application relates to a Hartmann wavefront sensor time super-resolution reconstruction method based on pixel-level three-spline interpolation with fusion of physical constraints and belongs to the technical field of beam wavefront measurement. BACKGROUND

[0002] As an important wavefront detection tool, the Shack-Hartmann Wavefront Sensor (SH-WFS) divides the incident wavefront through a microlens array and detects the centroid shift of each sub-aperture focal spot to reconstruct the wavefront phase. This technology is widely used in optical processing detection, adaptive optical systems, and laser beam quality analysis due to its non-contact, real-time, high precision, and other characteristics.

[0003] In actual dynamic measurement process, SH-WFS faces the inherent contradiction between time resolution and measurement signal-to-noise ratio. High-precision wavefront measurement usually requires a longer exposure time to ensure sufficient signal-to-noise ratio, but this will sacrifice the time resolution of the system, making it difficult to accurately capture the rapidly changing dynamic wavefront. Conversely, using short exposure time can improve the time resolution, but the lack of photons will lead to deterioration of the signal-to-noise ratio, seriously affecting the measurement accuracy.

[0004] The main shortcomings of existing time resolution improvement techniques include: 1. High-speed image sensor: high cost, limited by readout circuit speed and data transmission bandwidth; 2. Multi-sensor parallel sampling: complex system, requiring accurate time synchronization and complex calibration process; 3. Traditional linear interpolation method: simple calculation but limited accuracy, unable to accurately reconstruct complex wavefront dynamic characteristics.

[0005] These methods have their own shortcomings: hardware solutions are costly and complex to implement; linear interpolation methods can produce significant errors and cannot meet the needs of high-precision dynamic measurement. Therefore, an effective method is needed to achieve high-precision time super-resolution reconstruction under the condition of regular frame rate sampling. SUMMARY

[0006] In view of the deficiencies of the prior art, the present application provides a Hartmann wavefront sensor time super-resolution reconstruction method based on pixel-level three-spline interpolation fusing physical constraints, which detects and eliminates invalid points in the time sequence of each pixel in the preprocessing stage, eliminates non-physical mutations caused by noise, sampling anomalies or pixel clamping, and ensures the continuity and reliability of the input data. On this basis, a time three-spline interpolation model is constructed at the pixel level, and multiple levels of physical constraints are introduced in the interpolation process, so that the reconstructed high-time-resolution light spot sequence is not only mathematically continuous and smooth, but also physically more in line with the actual wavefront evolution characteristics. The present application can fully utilize the spatio-temporal correlation of the original light spot image, and realize higher precision and higher robustness of the wavefront time super-resolution reconstruction.

[0007] The present application adopts the following technical solutions: A Hartmann wavefront sensor time super-resolution reconstruction method based on pixel-level three-spline interpolation fusing physical constraints, comprising the following steps: (1) Single-CCD sequence input: receiving low-frame-rate light spot sequence data collected by a Hartmann wavefront sensor CCD; (2) Data preprocessing: for each pixel position, the intensity value sequence in the time dimension is extracted, the intensity value sequence is subjected to time domain median filtering for noise filtering, brightness normalization and time sequence alignment, and then invalid data points are detected and eliminated; (3) Three-spline interpolation modeling: a cubic spline interpolation model is independently constructed for the time sequence of each pixel after preprocessing, and physical constraints are introduced in the modeling process, including intensity range and energy consistency constraint, time change rate constraint and time smoothness and continuity constraint; (4) Coefficient solving: combining the physical constraint conditions, the coefficients of the three-spline interpolation function of each pixel are obtained by solving the improved three-diagonal equation set, so that the interpolation function meets the data fitting requirements and the optical physics law; (5) High-resolution sequence generation: new time sampling points are inserted between the original sampling time points according to the set super-resolution ratio, and the intensity value of each pixel at the new time point is calculated using the physical constraint three-spline interpolation function, to generate a high-time-resolution pixel-level time sequence; (6) Frame recombination: the interpolated pixel time sequence is recombined into a high-frame-rate image sequence.

[0008] Preferably, the specific implementation of step (2) comprises: (2.1) Time domain median filtering is used to remove impulse noise and maintain signal edge characteristics. First, the filtering window is a sliding time window, and the filtering window size is set to 3-7 frames, which can be adjusted according to the actual noise level and time resolution requirements.

[0009] Window sizeW satisfy: ; in, Represents the set of positive integers. k This is the half-width of the window, and its value is determined by the noise level. In temporal median filtering, for each pixel position Extract the intensity value sequence within the window along the time dimension, sort the intensity value sequence, and then take the median: ; in, Indicates the first in the image i line, number j The pixels in the column, at a given time point The intensity value after time-domain median filtering at each moment; This indicates the median operation, which involves sorting the sequence within the curly braces and then selecting the element at the middle position. They represent the first i line, number j The pixels in the column, at a given time point Raw intensity data at any given moment; Indicates the first line, number Column pixels in arrive After sorting the intensity values ​​at time 1, the one located at the 1st time... Elements of position; Edge processing is performed. For boundary frames at the beginning and end of the sequence, mirror expansion or truncation windows are used to ensure that all frames are effectively filtered. (2.2) Brightness fluctuations caused by exposure time variations are eliminated through brightness normalization. The frame with the most suitable exposure and the clearest light spot in the sequence is used as the reference frame. Usually, the middle frame of the sequence is selected or automatically selected through image quality evaluation. The brightness adjustment factor is calculated based on the statistical characteristics of the whole frame or local area. ; in, Indicates the first t The luminance scaling factor of the frame is used to scale the luminance of the first frame. t The brightness of the frame is adjusted to match the level of the reference frame; This represents the global average brightness of the reference frame, and its value is equal to... ; Indicates the first t The global average brightness of the frame is equal to ; M Indicates the number of pixels in the vertical direction; N Indicates the number of pixels in the horizontal direction; Indicates the first reference frame i line, number j The intensity value of the column pixel; Indicates the first t The first frame i line, number j Intensity values ​​of column pixels; This indicates a global averaging operation, which means averaging the intensity values ​​of all pixels in the image. For sequences with uneven brightness distribution, a block normalization strategy is adopted: ; in, Indicates the first t The first frame The brightness scaling factor for each block is used to adjust the brightness of that block to the same level as the corresponding block in the reference frame. Indicates the first reference frame The average brightness of each block is equal to the value of ; Indicates the first t The first frame The average brightness of each block is equal to the value of ; This represents the number of pixels in the horizontal direction. This represents the number of pixels in the vertical direction. Indicates the first The row pixel range of each block in the image; Indicates the first The column pixel range of each block in the image; Indicates the first t In frame The intensity value of the pixel at that location; This represents the average operation within the block, representing the operation on the th... The intensity values ​​of all pixels within the block are averaged. For cases with a large dynamic range, a nonlinear normalization method based on histogram matching is adopted to maintain the relative relationship of brightness distribution. (2.3) Ensure precise synchronization of timestamps for each frame through timing alignment, providing an accurate time reference for interpolation; extract precise exposure start and end times from the metadata of the camera acquisition system; for cases with acquisition delay, calculate the precise center time of each frame through timestamp interpolation or extrapolation: ; in, The timestamp representing the exposure center of a single frame image is the precise time reference for that frame in the time series, used for subsequent time alignment, interpolation, and other processing. This indicates the start time of the exposure process for that frame, recorded by the control system of the image acquisition device; represents the end time of the frame exposure process, also recorded by the control system of the acquisition device; represents the midpoint of the time interval from the beginning to the end of exposure, that is, the theoretical center time of the exposure process; represents the hardware and software delay of the acquisition system, used to correct the deviation between the theoretical center time and the actual effective exposure center, and improve the accuracy of the timestamp; Calculate the variance of the time interval of adjacent frames to detect whether there is frame loss or frame rate fluctuation: ; wherein, represents the variance of the time interval of adjacent frames, used to quantify the fluctuation degree of the time interval, and the larger the value, the worse the frame rate stability; N represents the total number of frames of the acquired image sequence; represents the time interval of the nth group of adjacent frames, that is, the timestamp difference value between the nth+1 frame and the nth frame; represents the average value of all adjacent frame time intervals, and the calculation formula is ; represents the timestamp of the nth+1 frame; represents the timestamp of the nth frame; The threshold is set in advance, and the threshold is usually 10%-20% of If is greater than the threshold, it is determined that there is frame rate fluctuation; if >2 , it is determined that there is frame loss between the adjacent frames; for frame rate fluctuation, the timestamp of the fluctuation frame is linearly interpolated and corrected, taking the timestamp of the adjacent normal frame as the reference to fit the reasonable timestamp of the fluctuation frame; for frame loss, according to the wavefront change rule before and after the frame loss, the pixel intensity value corresponding to the lost frame is completed by using cubic spline interpolation, so that the time continuity of the sequence is restored;

[0010] All timestamps of the frames are unified to a relative time coordinate system, taking the first frame as the time origin: ; represents the relative timestamp of the nth frame, that is, the time value after conversion taking the first frame as the time origin; represents the original timestamp of the nth frame, which is the original time data recorded by the acquisition system; represents the original timestamp of the first frame in the sequence, which is taken as the origin of the relative time coordinate system; N is the total number of frames of the acquired image sequence, and the value of n covers all frames to ensure that the timestamps of all frames are converted relative; (2.4) Assuming that the time sequence is: , each For each point, corresponding to an intensity value at a time point Calculate the mean value of the two points before and after: ; Then, calculate the difference between the point and the mean value : ; Use the standard deviation to measure the normal fluctuation range: ; Where is the mean value of the time series ; According to the standard deviation , set a threshold k , if is greater than , it is considered that is an invalid point; if it is judged that the point is invalid, replace its value with the mean value of the adjacent points before and after : .

[0011] Preferably, the implementation process of step (3) is: (3.1) Take the original sampling time point as the three spline nodes, and use the natural boundary condition, i.e. the second derivative of the two end nodes is zero: ; Where, represents the second derivative of the cubic spline interpolation function, representing the curvature change rate of the function at a certain time point; is the starting time point of the interpolation time node; is the end time point of the interpolation time node; At the beginning and end of the sequence, the wave front change rate is zero, which conforms to the actual physical scene of the starting and ending state; (3.2) For each pixel, construct a cubic polynomial on the time interval , and the standard three-spline form is: ; Where, represents the cubic spline interpolation function on the k th time sub-interval, corresponding to the pixel intensity interpolation model in the k th time sub-interval of the image sequence; t represents time, which is the input variable of the interpolation function; represents the starting time node of the k th time sub-interval, corresponding to the pixel intensity interpolation model in thek The timestamp of the frame is one of the nodes of the piecewise interpolation; respectively represent the coefficients of the cubic polynomial on the k th sub-interval, in order of the constant term, the linear term, the quadratic term, and the cubic term; k The index of the time sub-interval, taking values 0, 1,..., n-2 (n is the total number of time nodes, corresponding to the number of frames of the image sequence), represents that the n time nodes are divided into n-1 sub-intervals, each of which corresponds to an independent piecewise cubic polynomial; (3.3) Ensure the continuity of the function value, the first derivative, and the second derivative at the nodes; The continuity of the function value: ; The continuity of the first derivative: ; The continuity of the second derivative: ; The above continuity conditions constitute a linear equation system containing equations.

[0012] Preferably, the coefficient solving in step (4) adopts a tridiagonal matrix algorithm, the process of which is as follows: (4.1) According to the continuity conditions, construct a tridiagonal equation system about the second derivative : Wherein, ; denotes the length of the k th time sub-interval; denotes the first derivative of the cubic spline interpolation function at the th time node ; , , respectively denote the intensity value of the th pixel at the time node , which is the preprocessed pixel intensity sequence data; respectively denote the time nodes of the cubic spline interpolation, corresponding to the timestamps of different frames in the image sequence; denotes the first derivative of the cubic spline interpolation function of the th sub-interval at the time node ; (4.2) Eliminate the lower triangular elements by recursion, the goal being to transform the equation system containing lower triangular terms into an upper triangular form without lower triangular elements; For the second derivative the three-diagonal equation system is initialized, and let , the equation is directly simplified, and the coefficients after elimination are initialized: ; wherein, is the main diagonal, is the upper triangular, is the right end item; then the lower triangular items of the current equation are eliminated by using the simplification result of the step , and the elimination coefficients are updated: ; repeat this step until , and complete the elimination of the lower triangular elements; solve all node second-order derivatives by backward substitution: combine the boundary conditions , and recursively solve the second-order derivative of the last intermediate node from the last node: ; recursively solve from to : ; finally, the second-order derivatives of all nodes are obtained: ; (4.3) First, the intensity distribution characteristics of the pixel time series combined with physical constraints before interpolation are counted, including minimum value, maximum value, mean value and standard deviation; after interpolation, the newly generated pixel value is limited within the physical reasonable interval: ; wherein, is the standard deviation of the original pixel intensity sequence, representing the normal fluctuation amplitude of the intensity; at the same time, combined with the dynamic range of the detector, the interpolation result is further constrained to: ; avoid the appearance of physically impossible intensity values; In order to maintain the energy consistency of the light spot in the propagation process, the total light intensity of the interpolation frame is further constrained to change between adjacent frames, so as to maintain the physical reasonableness of the energy distribution, and then the time variation rate is constrained, the wavefront changes with limited bandwidth over time, and will not produce sudden intensity jump in adjacent sampling time; based on the original sampling sequence, the maximum time variation rate of each pixel is calculated: ; the change rate constraint is applied to the newly generated pixel value in the interpolation stage: ; ​This ensures that the interpolation result maintains the same rate of change as the original sequence, avoiding abrupt changes in intensity that do not conform to the laws of physical evolution; Combining physical constraints, for each time sub-interval k (k=0,1,...,n-2), using Calculate the coefficients: .

[0013] Preferably, the implementation process of step (5) is as follows: (5.1) Insert new time points uniformly between the original frames, with the insertion density determined by the super-resolution factor; let the target super-resolution factor be... In each original time interval Insertion A new point in time, a new point in time: ; in, Corresponding to the The first insertion within the original interval A new point in time; It is the first The length of each original time interval; the original N The frame sequence is obtained after super-resolution. frame; (5.2) For each new time point, calculate the pixel intensity value using the corresponding trispline function: ; (5.3) Ensure that the interpolation results are within a reasonable physical range. To avoid overflow exceptions; ; in, Indicates the pixel after clamping New time point The intensity value at that location; This represents the pixel obtained through interpolation, without physical constraint correction. New time point The original intensity value at the location; Indicates the first The first insertion within the original time interval A new point in time; The minimum and maximum pixel values ​​are determined by combining detector characteristics and optical system parameters to determine the reasonable range of pixel intensity. For the points needing clamping, the smoothing correction based on adjacent points is adopted, and for the clamped points, the intensity value of the point is corrected by using time distance weighted average with the intensity of adjacent original time points and the corrected intensity of adjacent new time points as references, so that the intensity of the corrected value is continuous with the intensity change of adjacent points, and sudden change is avoided.

[0014] The details of the present application can be seen from the prior art.

[0015] The beneficial effects of the present application are: In the preprocessing stage, the time sequence of each pixel is detected and removed for invalid points, and the non-physical mutation caused by noise, sampling anomaly or pixel clamping is eliminated to ensure the continuity and reliability of the input data. On this basis, a time three-spline interpolation model is constructed at the pixel level, and a multi-level physical constraint is introduced in the interpolation process, so that the reconstructed high-time-resolution light spot sequence is not only mathematically continuous and smooth, but also physically more in line with the actual wavefront evolution characteristics. Through the above innovative design, the present application can fully utilize the spatio-temporal correlation of the original light spot image, realize higher precision and higher robustness of the wavefront time super-resolution reconstruction, and has the following advantages:

[0016] Significant improvement in accuracy: the pixel-level three-spline interpolation fully utilizes the rich information of the original data, and the reconstruction accuracy is higher than that of linear interpolation; Maintain physical characteristics: the second-order smoothness of the three-spline ensures that the reconstructed wavefront sequence conforms to the physical laws of light wave propagation; Details are retained intact: the subtle changes in the shape of the light spot can be accurately reconstructed, and the edge definition is maintained; Good numerical stability: the three-spline interpolation algorithm is stable, and has certain inhibition ability to measurement noise; Optimization of computational efficiency: the three-diagonal matrix algorithm is adopted, and the computational complexity is low, which is suitable for real-time processing; Strong adaptability: suitable for various illumination conditions and noise level measurement environments. BRIEF DESCRIPTION OF DRAWINGS

[0017] The drawings accompanying the specification of this application form a part of this application and serve to further understand the application, the illustrative embodiments of the application and their descriptions serve to explain the application, and do not constitute an improper limitation on the application.

[0018] Figure 1 The flow chart of the Hartmann wavefront sensor time super-resolution reconstruction method based on pixel-level three-spline interpolation with fusion of physical constraints of the present application; Figure 2 The schematic diagram of the single-CCD low-frame-rate original light spot sequence; Figure 3 The schematic diagram of the high-frame-rate light spot sequence after three-spline interpolation reconstruction. DETAILED DESCRIPTION

[0019] In order to better understand the technical solutions in the specification, the technical solutions in the embodiments of the specification will be described clearly and completely below in combination with the drawings in the embodiments of the specification, but not limited to this. The unexplained embodiments of the present application are in accordance with the conventional techniques in the art.

[0020] Embodiment 1 A Hartmann wavefront sensor time super-resolution reconstruction method based on pixel-level three-spline interpolation combined with physical constraints, comprising the following steps: (1) Single-CCD sequence input: receiving low-frame-rate spot sequence data collected by a Hartmann wavefront sensor CCD; (2) Data preprocessing: for each pixel position, extract its intensity value sequence in the time dimension, perform noise filtering, brightness normalization and time sequence alignment processing on the intensity value sequence using time domain median filtering, and then detect and remove invalid data points; (3) Three-spline interpolation modeling: independently construct a cubic spline interpolation model for the time sequence of each pixel after preprocessing, and introduce physical constraints in the modeling process, including intensity range and energy consistency constraint, time change rate constraint, and time smoothness and continuity constraint; (4) Coefficient solving: combined with the physical constraint conditions, the coefficients of the three-spline interpolation function of each pixel are obtained by solving the improved three-diagonal equation set, so that the interpolation function meets the data fitting requirements and the optical physics law; (5) High-resolution sequence generation: new time sampling points are inserted between the original sampling time points according to the set super-resolution ratio, and the intensity value of each pixel at the new time point is calculated using the physical constraint three-spline interpolation function, to generate a high-time-resolution pixel-level time sequence; (6) Frame recombination: recombine the interpolated pixel time sequence into a high-frame-rate image sequence.

[0021] The specific implementation of step (2) includes: (2.1) Time domain median filtering is used to remove impulse noise and maintain signal edge characteristics. First, the filter window is a sliding time window, and the filter window size is set to 3-7 frames, which can be adjusted according to the actual noise level and time resolution requirements.

[0022] Window size W Satisfies: ; Wherein, Indicates the set of positive integers, k The window half-width is determined by the noise level; In the time domain median filtering, for each pixel position , where the intensity value sequence within the window is extracted in the time dimension, and the median value is taken after sorting the intensity value sequence: ; wherein, represents the intensity value of the pixel in the i th row and the j th column in the image after temporal median filtering at the time point ; represents the median operation, which refers to selecting the element in the middle position after sorting the sequence in the curly brackets; respectively represent the original intensity data of the pixel in the i th row and the j th column at the time point ; represents the element in the th position after sorting the intensity values of the pixel in the th row and the th column from the time point to the time point ; The edges are processed, and for the boundary frames at the beginning and end of the sequence, mirror extension or truncated window processing is adopted to ensure that all frames are effectively filtered; (2.2) Eliminate the brightness fluctuation caused by the change of exposure time by brightness normalization, take the frame with the most appropriate exposure and the clearest spot in the sequence as the reference frame, usually select the middle frame in the sequence or automatically select through image quality evaluation, calculate the brightness adjustment factor based on the statistical characteristics of the whole frame or local area: ; wherein, represents the brightness scaling coefficient of the t th frame, which is used to adjust the brightness of the t th frame to the same level as the reference frame; represents the global average brightness of the reference frame, whose value is equal to ; represents the global average brightness of the t th frame, whose value is equal to ; M represents the number of pixels in the vertical direction; N represents the number of pixels in the horizontal direction; represents the intensity value of the pixel in the i th row and the j th column in the reference frame; represents the intensity value of the pixel in the t th row and the i th column in the j th frame; represents global average operation, which means taking average of intensity values of all pixels in the image; For sequences with uneven brightness distribution, block normalization strategy is adopted: ; where, represents the brightness scaling factor of the t th block in the th frame, which is used to adjust the brightness of the block to the same level as the corresponding block in the reference frame; represents the average brightness of the th block in the reference frame, which is equal to ; represents the average brightness of the t th block in the th frame, which is equal to ; is the number of pixels in horizontal direction; is the number of pixels in vertical direction; represents the row pixel range of the th block in the image; represents the column pixel range of the th block in the image; represents the intensity value of the pixel at the t th position in the th frame; represents block-wise average operation, which means taking average of intensity values of all pixels in the th block; For cases with large dynamic range, histogram matching based nonlinear normalization method is adopted to maintain the relative relationship of brightness distribution; (2.3) By time series alignment, ensure the accuracy of timestamp synchronization of each frame, provide accurate time reference for interpolation; Extract the accurate exposure start time and end time from the metadata of the camera acquisition system, for cases with acquisition delay, calculate the accurate center time of each frame by timestamp interpolation or extrapolation: ; where, represents the exposure center timestamp of a single frame image, which is the accurate time reference of the frame in time series, used for subsequent time series alignment, interpolation and other processing; represents the start time of the exposure process of the frame, recorded by the control system of the image acquisition device; represents the end time of the exposure process of the frame, also recorded by the control system of the acquisition device; represents the midpoint of the time interval from the start to the end of exposure, that is, the theoretical center time of the exposure process; Hardware and software delay of the acquisition system, used to correct the deviation between the theoretical center time and the actual effective exposure center, and improve the accuracy of the timestamp; Calculate the variance of the adjacent frame time interval to detect whether there is frame loss or frame rate fluctuation: ; Wherein, The variance of the adjacent frame time interval is used to quantify the fluctuation degree of the time interval, and the larger the value is, the worse the frame rate stability is; N Total frame number of the acquired image sequence; The time interval of the nth group of adjacent frames, that is, the timestamp difference value between the nth+1 frame and the nth frame; The average value of all adjacent frame time intervals, the calculation formula is ; The timestamp of the nth+1 frame; The timestamp of the nth frame; The threshold is set in advance, and the threshold is usually 10%-20% of If is greater than the threshold, it is determined that there is frame rate fluctuation; if >2 , it is determined that there is frame loss between the adjacent frames; for frame rate fluctuation, the timestamps of the fluctuation frames are linearly interpolated and corrected, taking the timestamps of the adjacent normal frames as the reference to fit the reasonable timestamps of the fluctuation frames; for frame loss, according to the wavefront change rule before and after the frame loss, the pixel intensity value corresponding to the lost frame is completed by using cubic spline interpolation, so that the time continuity of the sequence is restored;

[0023] Unify the timestamps of all frames to a relative time coordinate system, taking the first frame as the time origin: ; The relative timestamp of the nth frame, that is, the time value after conversion taking the first frame as the time origin; The original timestamp of the nth frame, which is the original time data recorded by the acquisition system; The original timestamp of the first frame in the sequence, which is taken as the origin of the relative time coordinate system; N Total frame number of the acquired image sequence, the value range of n covers all frames, ensuring that the timestamps of all frames are converted relative; (2.4) Assuming that the time sequence is: , each corresponds to an intensity value at a time point, and for each point Calculate the average value of the two points before and after: ; Then, calculate the point The difference from the mean value : ; The normal fluctuation range is measured by the standard deviation ; Wherein is the mean value of the time series ; A threshold is set according to the standard deviation k If is greater than , it is considered that is an invalid point; if it is judged that the point is invalid, its value is replaced by the mean value of the adjacent points before and after : .

[0024] The implementation process of step (3) is as follows: (3.1) The original sampling time point is taken as the three-spline node, and the natural boundary condition is adopted, that is, the second derivative of the two end nodes is zero: ; Wherein, represents the second derivative of the cubic spline interpolation function, representing the curvature change rate of the function at a certain time point; is the starting time point of the interpolation time node; is the terminal time point of the interpolation time node; At the beginning and end of the sequence, the wave front change rate is zero, which conforms to the starting and ending state of the actual physical scene; (3.2) For each pixel, a cubic polynomial is constructed on the time interval The standard three-spline form is as follows: ; Wherein, represents the cubic spline interpolation function on the k th time sub-interval, corresponding to the pixel intensity interpolation model in the k th time sub-interval in the image sequence; t represents the time, which is the input variable of the interpolation function; represents the starting time node of the k th time sub-interval, corresponding to the time stamp of the k th frame in the image sequence, which is one of the nodes of the segmented interpolation; respectively represent the coefficients of the cubic polynomial on the k th sub-interval, which are the constant term, the first-order term, the second-order term and the third-order term coefficients, respectively;​k index of time sub-interval, value range 0, 1, …, n-2 (n is the total number of time nodes, corresponding to the frame number of image sequence), indicating that n time nodes are divided into n-1 sub-intervals, each sub-interval corresponds to an independent cubic polynomial segment; (3.3) Ensure the continuity of function value, first derivative and second derivative at the node; Function value continuity: ; First derivative continuity: ; Second derivative continuity: ; The above continuity conditions constitute a linear equation system containing equations.

[0025] The coefficient solving in step (4) adopts a three-diagonal matrix algorithm, the process is as follows: (4.1) According to the continuity condition, a three-diagonal equation group about the second derivative is constructed: ; Wherein, ; denotes the length of the k time sub-interval; denotes the first derivative of the cubic spline interpolation function at the time node ; , , respectively denote the intensity value of the pixel at the time node , which is the preprocessed pixel intensity sequence data; respectively denote the time nodes of cubic spline interpolation, corresponding to the time stamp of different frames in the image sequence; denotes the first derivative of the cubic spline interpolation function of the sub-interval at the time node ; (4.2) Eliminate the lower triangular elements by recursion, the goal is to transform the equation group containing lower triangular elements into an upper triangular form without lower triangular elements; For the three-diagonal equation group about the second derivative in step (4.1), initialization is carried out, let , , directly simplify the equation, and the initialized eliminated coefficient is: ; where, is the main diagonal, is the upper triangle, is the right end item; then use the simplified result of step to eliminate the lower triangle item of the current equation , update the elimination coefficient: ; repeat this step until , complete the elimination of the lower triangle element; solve all node second-order derivatives by backward substitution Combine the boundary conditions , recursively from the last node to the last intermediate node to solve the second-order derivative of the last intermediate node: ; From , recursively to : ; Finally, the second-order derivative of all nodes is obtained: ; (4.3) First, the intensity distribution characteristics of the pixel time series combined with physical constraints before interpolation are calculated, including minimum value, maximum value, mean value and standard deviation; after interpolation, the newly generated pixel value is limited within the physical reasonable interval: ; where, is the standard deviation of the original pixel intensity sequence, representing the normal fluctuation amplitude of the intensity; at the same time, combined with the dynamic range of the detector, the interpolation result is further constrained to: ; avoid the appearance of physically impossible intensity values; In order to maintain the energy consistency of the light spot in the propagation process, further constrain the total light intensity of the interpolation frame to change between adjacent frames within a limited range, so as to maintain the physical reasonableness of the energy distribution, and then constrain the time variation rate, the wavefront changes with time has a limited bandwidth, and will not produce sudden intensity jump in adjacent sampling time; based on the original sampling sequence, calculate the maximum time variation rate of each pixel: ; In the interpolation stage, the change rate constraint is applied to the newly generated pixel value: ; so as to ensure that the interpolation result maintains the same change speed as the original sequence, avoiding intensity mutation that does not conform to the physical evolution law; Combined with physical constraints, for each time sub-interval k (k=0,1,...,n-2), use Calculate the coefficients: .

[0026] The implementation process of step (5) is as follows: (5.1) Insert new time points uniformly between the original frames, with the insertion density determined by the super-resolution factor; let the target super-resolution factor be... In each original time interval Insertion A new point in time, a new point in time: ; in, Corresponding to the The first insertion within the original interval A new point in time; It is the first The length of the original time interval; the original N The frame sequence is obtained after super-resolution. frame; (5.2) For each new time point, calculate the pixel intensity value using the corresponding trispline function: ; (5.3) Ensure that the interpolation results are within a reasonable physical range. To avoid overflow exceptions; ; in, Indicates the pixel after clamping New time point The intensity value at that location; This represents the pixel obtained through interpolation, without physical constraint correction. New time point The original intensity value at the location; Indicates the first The first insertion within the original time interval A new point in time; The minimum and maximum pixel values ​​are determined by combining detector characteristics and optical system parameters to determine the reasonable range of pixel intensity. For points that need to be clamped, a smooth correction based on neighboring points is adopted. For clamped points, the intensity value of the point is corrected by using the intensity of the adjacent original time point and the intensity of the corrected adjacent new time point as references, and the intensity value of the point is corrected by time distance weighted average, so that the corrected value is continuous with the intensity change of the adjacent points and avoids introducing abrupt changes.

[0027] Example 2 A temporal super-resolution reconstruction method for Hartmann wavefront sensors based on pixel-level trispline interpolation, incorporating physical constraints, is proposed. Figure 1As shown, first, a sequence of light spots is collected by a single CCD sensor at a certain frame rate (such as Figure 2 ), and after preprocessing such as time domain median filtering, brightness normalization and time sequence alignment, invalid point removal, a stable low frame rate sequence is obtained. The rapid dynamic change of the wavefront phase can be observed in the original sequence, but due to the frame rate limitation, the complete evolution process cannot be accurately captured, resulting in the loss of dynamic information. At this time, if the low frame rate sequence is directly used for wavefront reconstruction, the rapidly changing wavefront characteristics cannot be accurately reflected, resulting in a loss of dynamic measurement accuracy. After preprocessing, the sequence enters the pixel-level three-spline interpolation processing stage. A three-spline function model in the time dimension is independently constructed for each pixel. The interpolation coefficients are solved by a three-diagonal matrix algorithm combined with physical constraints. New time points are inserted between the original frames, and the pixel intensity values are calculated. Finally, a high time resolution sequence with a frame rate of k times is generated (such as Figure 3 ). From the processing results, it can be seen that the complete dynamic evolution process of the wavefront phase is effectively reconstructed by the method, solving the problem of dynamic information loss caused by insufficient sampling rate, and providing a reliable technical means for high-precision dynamic wavefront measurement.

[0028] (1) A conventional Hartmann wavefront sensor is used to collect a sequence of light spot data at a certain frame rate. The image resolution is 128x128 pixels. Continuous collection is used as the input data source to provide basic data support for subsequent time super-resolution processing.

[0029] (2) The original sequence collected is subjected to time domain median filtering processing, and a 3-frame window is used to effectively suppress impulse noise. The brightness fluctuation of the subsequent frames is eliminated by brightness normalization with the first frame as the reference. Time sequence alignment is performed based on the camera accurate timestamp to ensure the accuracy of the time reference, and then invalid points are removed to lay a foundation for subsequent interpolation processing.

[0030] For each pixel position (i,j) in the image, the intensity value change sequence in the entire time dimension is independently extracted, and a total of 128x128=16384 independent time sequences are obtained, each containing 100 time point intensity value data, providing complete input data for pixel-level time interpolation.

[0031] (3) A three-spline function model is constructed for each pixel's independent time sequence. The original sampling time points are used as interpolation nodes, natural boundary conditions are used to ensure that the second derivatives at the beginning and end of the sequence are zero, and strict continuity conditions are used to ensure smooth transition of the function at the nodes, establishing a complete mathematical interpolation framework.

[0032] (4) Combined with physical constraints, the high efficient tri-diagonal matrix algorithm is used to solve the coefficients, first, the tri-diagonal equations about the second derivative are constructed, the intermediate variables are calculated through forward elimination, then the backward substitution is carried out to obtain the second derivative values of all nodes, finally, all coefficients of the tri-cubic function are calculated according to the values.

[0033] (5) Realize k time super-resolution reconstruction, a plurality of new time points are inserted in each original time interval, the pixel intensity values of the new time points are calculated by using the obtained tri-cubic function, and the original low frame sequence is expanded to k times high time resolution sequence.

[0034] (6) The interpolated pixel time sequence is recombined into a complete image sequence to form a three-dimensional data cube, and finally the frame rate is increased to a high k times high time resolution spot sequence, and the time super-resolution processing flow is completed.

[0035] The above describes the preferred embodiments of the present application, and it should be noted that for those skilled in the art, without departing from the principles of the present application, a number of improvements and refinements can be made, and these improvements and refinements should also be considered as the protection scope of the present application.

Claims

1. A Hartmann wavefront sensor temporal super-resolution reconstruction method based on pixel-level trispline interpolation and incorporating physical constraints, characterized in that, Includes the following steps: (1) Single CCD sequence input: Receives low frame rate spot sequence data acquired by Hartmann wavefront sensor CCD; (2) Data preprocessing: For each pixel location, extract its intensity value sequence in the time dimension, and perform noise filtering, brightness normalization and time alignment on the intensity value sequence by temporal median filtering. Then, invalid data points are detected and removed. (3) Trispline interpolation modeling: Construct a cubic spline interpolation model independently for the time series of each preprocessed pixel, and introduce physical constraints in the modeling process; (4) Coefficient solution: Combining physical constraints, the coefficients of the trispline interpolation function of each pixel are obtained by solving the improved tridiagonal equation system, so that the interpolation function can simultaneously meet the data fitting requirements and the laws of optical physics. (5) High-resolution sequence generation: Insert new time sampling points between the original sampling times according to the set super-resolution ratio, and use the physical constraint trispline interpolation function to calculate the intensity value of each pixel at the new time point to generate a pixel-level time series with high time resolution; (6) Frame reconstruction: Recombining the interpolated pixel time series into a high frame rate image sequence.

2. The Hartmann wavefront sensor temporal super-resolution reconstruction method based on pixel-level trispline interpolation with fused physical constraints as described in claim 1, characterized in that, The specific implementation of step (2) includes: (2.1) Impulse noise is removed by time-domain median filtering. First, a sliding time window is used for filtering, and the window size is... W satisfy: ； in, Represents the set of positive integers. k The width is half the width of the window; In temporal median filtering, for each pixel position Extract the intensity value sequence within the window along the time dimension, sort the intensity value sequence, and then take the median: ； in, Indicates the first in the image i line, number j The pixels in the column, at a given time point The intensity value after time-domain median filtering at each moment; This indicates the median operation, which involves sorting the sequence within the curly braces and then selecting the element at the middle position. They represent the first i line, number j The pixels in the column, at a given time point Raw intensity data at any given moment; Indicates the first line, number Column pixels in arrive After sorting the intensity values ​​at time 1, the one located at the 1st time... Elements of position; Edge processing is performed. For boundary frames at the beginning and end of the sequence, mirror expansion or truncation windows are used to ensure that all frames are effectively filtered. (2.2) Brightness fluctuations caused by exposure time variations are eliminated through brightness normalization. The frame with the most suitable exposure and the clearest light spot in the sequence is used as the reference frame. The brightness adjustment factor is calculated based on the statistical characteristics of the entire frame or a local area. ； in, Indicates the first t The luminance scaling factor of the frame is used to scale the luminance of the first frame. t The brightness of the frame is adjusted to match the level of the reference frame; This represents the global average brightness of the reference frame, and its value is equal to... ; Indicates the first t The global average brightness of the frame is equal to ; M Indicates the number of pixels in the vertical direction; N Indicates the number of pixels in the horizontal direction; Indicates the first reference frame i line, number j The intensity value of the column pixel; Indicates the first t The first frame i line, number j Intensity values ​​of column pixels; This indicates a global averaging operation, which means averaging the intensity values ​​of all pixels in the image. For sequences with uneven brightness distribution, a block normalization strategy is adopted: ； in, Indicates the first t The first frame The brightness scaling factor for each block is used to adjust the brightness of that block to the same level as the corresponding block in the reference frame. Indicates the first reference frame The average brightness of each block is equal to the value of ; Indicates the first t The first frame The average brightness of each block is equal to the value of ; This represents the number of pixels in the horizontal direction. This represents the number of pixels in the vertical direction. Indicates the first The row pixel range of each block in the image; Indicates the first The column pixel range of each block in the image; Indicates the first t In frame The intensity value of the pixel at that location; This represents the average operation within the block, representing the operation on the th... The intensity values ​​of all pixels within the block are averaged. For cases with a large dynamic range, a nonlinear normalization method based on histogram matching is adopted to maintain the relative relationship of brightness distribution. (2.3) Ensure precise synchronization of timestamps for each frame through timing alignment, providing an accurate time reference for interpolation; extract precise exposure start and end times from the metadata of the camera acquisition system; for cases with acquisition delay, calculate the precise center time of each frame through timestamp interpolation or extrapolation: ； in, The timestamp representing the exposure center of a single frame image is the precise time reference for that frame in the time series, used for subsequent time alignment, interpolation, and other processing. This indicates the start time of the exposure process for that frame, recorded by the control system of the image acquisition device; This indicates the end time of the exposure process for that frame, which is also recorded by the control system of the acquisition device; This indicates the hardware and software latency of the data acquisition system; Calculate the variance of the time interval between adjacent frames to detect whether there are dropped frames or frame rate fluctuations: ； in, It represents the variance of the time interval between adjacent frames, used to quantify the degree of fluctuation in the time interval. The larger the value, the worse the frame rate stability. N This indicates the total number of frames in the acquired image sequence; This represents the time interval between the nth group of adjacent frames, that is, the difference in timestamps between the (n+1)th and nth frames; This represents the average time interval between all adjacent frames; This represents the timestamp of the (n+1)th frame; This represents the timestamp of the nth frame; Pre-set threshold, if If the value exceeds the threshold, a frame rate fluctuation is determined to exist; if >2 Then it is determined that there is a frame drop between adjacent frames in the group; for frame rate fluctuations, the timestamp of the fluctuating frame is corrected by linear interpolation, and the reasonable timestamp of the fluctuating frame is fitted with the timestamp of the adjacent normal frame as the reference; for frame drops, according to the wavefront change law of the frame before and after the frame drop, cubic spline interpolation is used to complete the pixel intensity value corresponding to the frame drop, so as to restore the temporal continuity of the sequence. Unify the timestamps of all frames to a relative time coordinate system, with the first frame as the time origin: ； This represents the relative timestamp of the nth frame, which is the time value after conversion with the first frame as the time origin; This represents the original timestamp of the nth frame, which is the original time data recorded by the acquisition system; Represents the original timestamp of the first frame in the sequence, serving as the origin of the relative time coordinate system; N The total number of frames in the acquired image sequence is n, and the range of n values ​​covers all frames to ensure that the timestamps of all frames are converted to relative values. (2.4) Assume the time series is as follows: Each For each point in time, the intensity value corresponds to that point. Calculate the mean of the two points before and after: ； Then, calculate the points. with the mean Differences: ； Use standard deviation to measure the normal range of fluctuation: ； in The mean of the time series is... ; According to standard deviation Set a threshold k ,if Greater than Then it is believed It is an invalid point; if the point is determined... It is invalid; replace its value with the mean of its preceding and following neighbors. : 。 3. The Hartmann wavefront sensor temporal super-resolution reconstruction method based on pixel-level trispline interpolation with fused physical constraints as described in claim 2, characterized in that... In step (2.1), the filter window size is set to 3-7 frames.

4. The Hartmann wavefront sensor temporal super-resolution reconstruction method based on pixel-level trispline interpolation with fused physical constraints as described in claim 2, characterized in that, In step (2.3), the threshold is taken as follows: 10% to 20%.

5. The Hartmann wavefront sensor temporal super-resolution reconstruction method based on pixel-level trispline interpolation with fused physical constraints as described in claim 1, characterized in that, The implementation process of step (3) is as follows: (3.1) Using the original sampling time points as the nodes of the trispline, and applying natural boundary conditions, i.e., the second derivatives of the two endpoints are zero: ； in, The second derivative of the cubic spline interpolation function represents the rate of change of the curvature of the function at a certain time point; This is the starting time point for the interpolation time node; This is the end time point of the interpolation time node; At the beginning and end of the sequence, the wavefront change rate is zero, which is consistent with the start and end states of the actual physical scenario. (3.2) For each pixel, in the time interval Constructing a cubic polynomial, the standard trispline form is: ； in, Indicates the first k The cubic spline interpolation function over the nth time sub-interval corresponds to the nth time interval in the image sequence. k Pixel intensity interpolation model within a time sub-interval; t Representing time, it is the input variable of the interpolation function; Indicates the first k The starting time node of the nth time sub-interval corresponds to the nth time in the image sequence. k The timestamp of a frame is one of the nodes in the segmented interpolation; They represent the first k The coefficients of the cubic polynomial over each subinterval are, in order, the coefficients of the constant term, the linear term, the quadratic term, and the cubic term; k The index of the time sub-interval, with a value range of 0, 1, ..., n-2, indicates that n time nodes are divided into n-1 sub-intervals, and each sub-interval corresponds to an independent piecewise cubic polynomial; (3.3) Ensure that the function values, first derivatives, and second derivatives are continuous at the nodes; Function values ​​are continuous: ； The first derivative is continuous: ； The second derivative is continuous: ； The above continuity conditions constitute the inclusion A system of linear equations.

6. The Hartmann wavefront sensor temporal super-resolution reconstruction method based on pixel-level trispline interpolation with fused physical constraints as described in claim 5, characterized in that, The coefficients in step (4) are solved using a tridiagonal matrix algorithm, as follows: (4.1) Construct the second derivative based on the continuity condition. The tridiagonal system of equations: ； in, ; Indicates the first k The length of each time sub-interval; This indicates that the cubic spline interpolation function is at the th... Time nodes The first derivative at that point; , , They represent the first Each pixel at the time node The intensity value at that location; These represent the time points of the cubic spline interpolation; Indicates the first The cubic spline interpolation function for each sub-interval at time node The first derivative at that point; (4.2) Eliminate the lower triangular elements by recursive calculation, with the goal of transforming the system of equations containing lower triangular terms into upper triangular form without lower triangular elements; Regarding the second derivative in step (4.1) The tridiagonal system of equations is initialized by letting... , The equation is simplified directly, and the coefficients after initial elimination are: ； in, Main diagonal It is an upper triangle. For the right-hand term; then use The simplification step eliminates the lower triangular terms of the current equation. Update the elimination coefficients: ； Repeat this step until Complete the elimination of lower triangular elements; solve for the second derivative of all nodes by back substitution. Combined with boundary conditions Starting from the last node and working backwards, we can solve for the second derivative of the last intermediate node: ； from Reverse recursion to : ； Finally, the second derivatives of all nodes are obtained: ; (4.3) First, the intensity distribution characteristics of the pixel time series before interpolation, combined with physical constraints, are statistically analyzed, including the minimum, maximum, mean, and standard deviation; after interpolation, the newly generated pixel values ​​are restricted to this physically reasonable range: ； in, The standard deviation of the original pixel intensity sequence represents the normal fluctuation range of intensity; simultaneously, combined with the detector's dynamic range, the interpolation results are further constrained to: ； To avoid physically impossible strength values; To maintain the energy consistency of the light spot during propagation, the total light intensity of the interpolated frames is further constrained to have limited variation between adjacent frames, thus preserving the physical rationality of the energy distribution. Then, the temporal rate of change is constrained, ensuring that the wavefront's change over time has a finite bandwidth, preventing sudden intensity jumps between adjacent sampling times. The maximum temporal rate of change for each pixel is calculated based on the original sampling sequence. ； Apply a rate-of-change constraint to the newly generated pixel values ​​during the interpolation stage: ； This ensures that the interpolation result maintains the same rate of change as the original sequence, avoiding abrupt changes in intensity that do not conform to the laws of physical evolution; Combining physical constraints, for each time sub-interval k (k=0,1,...,n-2), calculate the coefficients: 。 7. The Hartmann wavefront sensor temporal super-resolution reconstruction method based on pixel-level trispline interpolation with fused physical constraints as described in claim 6, characterized in that, The implementation process of step (5) is as follows: (5.1) Insert new time points uniformly between the original frames, with the insertion density determined by the super-resolution factor; let the target super-resolution factor be... In each original time interval Insertion A new point in time, a new point in time: ； in, Corresponding to the The first insertion within the original interval A new point in time; It is the first The length of the original time interval; the original N The frame sequence is obtained after super-resolution. frame; (5.2) For each new time point, calculate the pixel intensity value using the corresponding trispline function: ； (5.3) Ensure that the interpolation results are within a reasonable physical range. To avoid overflow exceptions; ； in, Indicates the pixel after clamping New time point The intensity value at that location; This represents the pixel obtained through interpolation, without physical constraint correction. New time point The original intensity value at the location; Indicates the first The first insertion within the original time interval A new point in time; These are the minimum and maximum pixel values ​​determined jointly by the detector characteristics and the optical system parameters; For points that need to be clamped, a smooth correction based on neighboring points is adopted. For clamped points, the intensity value of the point is corrected by using the intensity of the adjacent original time point and the intensity of the corrected adjacent new time point as references, and the intensity value of the point is corrected by time distance weighted average, so that the corrected value is continuous with the intensity change of the adjacent points and avoids introducing abrupt changes.

Images