Quasi-Gaussian pulse forming method based on cubic spline interpolation

Abstract

The invention belongs to the technical field of photon counting imaging, and particularly relates to a quasi-Gaussian pulse forming method based on cubic spline interpolation. Comprising the following steps: S1, performing pole-zero processing on an exponential decay signal to obtain a pole-zero cancellation signal; s2, fitting the pole zero cancellation signal to obtain a fitting parameter; solving the pole-zero cancellation signal based on a deconvolution recursion formula to obtain a deconvolution signal; s3, cubic spline interpolation fitting calculation is carried out on the standard Gaussian window, and an interpolation curve formula set of the Gaussian pulse shape is designed; s4, calculating the interpolation curve formula group to obtain a forming transfer function recursion formula; and S5, inputting the deconvolution signal into a forming transfer function recursion formula, and outputting a quasi-Gaussian forming signal. The method is high in universality, high in stability, adjustable in parameter and flexible in design.

Description

Technical Field

[0001] The invention belongs to the technical field of photon counting imaging, and in particular relates to a quasi-Gaussian pulse shaping method based on cubic spline interpolation. Background Art

[0002] In the field of imaging observation of weak light intensity targets in space, the photon counting imaging detector method based on microchannel plate (MCP) is widely used. The photon counting imaging detector mainly includes a photodetector based on a position-sensitive anode, a charge-sensitive amplifier circuit, a pulse shaper, a pulse peak identification and an image inversion module. For the exponential decay pulse sequence output by the preamplifier, the existing solution is to use a pulse shaper to shape the voltage pulse output by the preamplifier to change the shape of the pulse to facilitate subsequent signal acquisition. From the perspective of the frequency domain, the shaper can also be used as a filter to attenuate noise and improve the signal-to-noise ratio by limiting the bandwidth. The commonly used shaper is CR-RC n Shapers and Sallen-Key shapers, these two shapers can make the output pulse top have a more gentle asymmetric Gaussian-shaped pulse than the original input signal, and can achieve more accurate pulse amplitude extraction. And these shapers are all implemented using hardware analog circuits. In addition, there is another method in the prior art that performs pole-zero phase cancellation processing on the preamplifier signal, digitizes it through an analog-to-digital converter, and then uses a processor such as FPGA to implement trapezoidal shaping of the pulse signal, and finally extracts the amplitude of the trapezoidal shaped signal. Since the trapezoidal pulse has a smoother pulse top, more accurate pulse amplitude extraction can be achieved.

[0003] However, the existing technology has the following disadvantages: 1.CR-RC n The shaper and Sallen-Key shaper require a lot of hardware components to realize pulse shaping, and the width and shape of the shaped pulse cannot be adjusted flexibly and adaptively, and the peak extraction cannot achieve high accuracy. Although the shaper constructed by hardware analog circuit can realize pulse filtering and Gaussian-like shaping, due to the fixed circuit shaping parameters, the width and shape of the shaped pulse cannot be adjusted adaptively, and cannot meet the application requirements of high resolution and high counting rate.

[0004] 2. Compared with the Sallen-Key shaper, the trapezoidal shaping method can flexibly design the pulse width, and the time domain characteristics of the pulse are more conducive to the extraction of the pulse amplitude. However, due to the limitations of the trapezoidal design itself, the frequency domain filtering effect is not as good as the quasi-Gaussian filtering method with a smoother edge roll-off characteristic, which makes it impossible to better suppress high-frequency signals and the peak extraction cannot achieve a high accuracy. Summary of the invention

[0005] In view of this, the present invention aims to provide a quasi-Gaussian pulse shaping method based on cubic spline interpolation to solve the shortcomings of the prior art that the width and shape of the shaped pulse cannot be adaptively and flexibly adjusted, and the peak extraction cannot achieve high precision. The present invention combines the convenient time domain design characteristics of the digital shaping method and the high-frequency suppression characteristics of Gaussian shaping, and realizes adaptive adjustment of the pulse width and shape by setting different control points. The filtering effect of the present invention is closer to the ideal Gaussian window. The present invention has high versatility, strong stability, adjustable parameters and flexible design.

[0006] To achieve the above object, the technical solution created by the present invention is implemented as follows: A quasi-Gaussian pulse shaping method based on cubic spline interpolation specifically comprises the following steps: S1: Acquire an exponential decay signal, and perform pole-zero processing on the exponential decay signal to obtain a pole-zero cancellation signal; S2: fitting the pole-zero cancellation fragment using the three-exponential model to obtain fitting parameters; recursively deducing the deconvolution system transfer function based on the fitting parameters to obtain the deconvolution recursive formula, and solving the pole-zero cancellation signal based on the deconvolution recursive formula to obtain the deconvolution signal; S3: setting a standard Gaussian window, performing cubic spline interpolation fitting calculation on the standard Gaussian window, and designing an interpolation curve formula group of the Gaussian pulse shape; S4: extrapolating the interpolation curve formula group to obtain a recursive formula of the forming transfer function; S5: Input the deconvolution signal into the shaping transfer function recursive formula and output a quasi-Gaussian shaping signal.

[0007] Furthermore, step S1 specifically includes the following steps: S11: Input signal It is represented as an exponential decay signal, and the Laplace transform of the exponential decay signal is performed: ; ; in, is an input signal in exponential decay form, is the time interval of the exponential decay signal from the pulse start time to the time when the pulse decays to 1 / e of the peak value, is a natural constant, t is time, s is a complex variable, is the exponential decay signal after Laplace transform; S12: Set the first time decay constant , the second time decay constant , using the first time decay constant and the second time decay constant to calculate the pole-zero cancellation transfer function in the s domain : ; Where s is a complex variable; S13: Convert the pole-zero cancellation transfer function to the time domain and discretize it to obtain the pole-zero cancellation recursive expression. Discretize to obtain an exponential decay signal , the exponential decay signal Substitute the pole-zero cancellation recursive expression to obtain the pole-zero cancellation signal : ; ; ; ; ; in, is the sampling period, , b, c and d are all intermediate variables for calculating the pole-zero cancellation signal and have no physical meaning. n is the serial number of the discrete time point.

[0008] Furthermore, the three-index model is: ; in, It is a three-exponential model, A is the first fitting parameter, B is the second fitting parameter, C is the third fitting parameter, D is the fourth fitting parameter, E is the fifth fitting parameter, F is the sixth fitting parameter, A, B, C ≥ 0, A+B+C=0, |A|>|B|>|C|.

[0009] Furthermore, step S2 specifically includes the following steps: S21: Cancellation signal at pole zero Select the pole-zero cancellation fragment , the pole-zero cancellation fragments were fitted using a three-exponential model to obtain fitting parameters; S22: Substitute the fitting parameters into the three-exponential model used in step S2 for Z transformation and organize them into Rational fractions of S23: recursively deducing the deconvolution system transfer function based on the rational fraction in step S22 to obtain a deconvolution recursive formula; S24: Input the pole-zero cancellation signal into a deconvolution recursive formula for solving, and obtain a deconvolution signal.

[0010] Furthermore, the pole-zero cancellation segment is a complete, independent and non-pile-up pulse signal selected from the pole-zero cancellation signal, and the pulse peak amplitude of the pulse signal is at least 10 times the noise standard deviation.

[0011] Furthermore, the deconvolution system transfer function of the three-exponential model is: ; in, is a three-index model, is the impact function in the z domain; The deconvolution recursive formula of the three-exponential model is:

[0012] ; ; ; ; ; ; ; ; ; ; in, , , , , , , , and are the coefficients of the deconvolution system transfer function expressing the three-exponential model and have no physical meaning. is the deconvolution signal.

[0013] Furthermore, step S3 specifically includes the following steps: S31: Set the standard Gaussian window w: ; ; Among them, the length of the standard Gaussian window is L-1, L is the number of sampling points, L is an odd number, , is the width factor, is the standard deviation of the Gaussian curve, and n is the serial number of the discrete time point; S32: among the L sampling points, N sampling points are selected as interpolation control points, and the N sampling points are selected symmetrically with the sampling point n0=(L-1) / 2+1 as the center; S33: Fitting N interpolation control points based on the cubic sample interpolation method to obtain a fitting curve, and the N interpolation control points divide the fitting curve into N segments, and calculating the interpolation curve formula group of the N segments of the fitting curve: ; ;

[0014] ;

[0015] ; ; in, is the first fitting curve, is the second fitting curve, is the i-th segment fitting curve, is the N-1th segment fitting curve, is the Nth segment fitting curve, is the x-coordinate of the i-th point among the N control points, The polynomial coefficients of the i-th curve are obtained by cubic spline fitting using control points.

[0016] Furthermore, step S4 specifically includes the following steps: S41: discretize the interpolation curve formula group, and perform a neighboring term subtraction operation on the discretized interpolation curve formula group: ; ;

[0017] ;

[0018] ; ; in, is the first term of Gaussian pulse shaping fitted by cubic spline interpolation, The second term of Gaussian pulse shaping fitted by cubic spline interpolation, is the i-th term of Gaussian pulse shaping with cubic spline interpolation, is the N-1th term of Gaussian pulse shaping with cubic spline interpolation, is the Nth term of Gaussian pulse shaping with cubic spline interpolation, is the first segment of the cubic spline fitting curve after discretization, is the second segment of the cubic spline fitting curve after discretization, is the i-1th segment cubic spline fitting curve after discretization, is the discretized cubic spline fitting curve of the ith segment, is the N-2th segment of the discretized cubic spline fitting curve, is the N-1th segment of the discretized cubic spline fitting curve, For A step function where the position changes stepwise; S42: Use long division to To sort out and sort out Perform a Z transform: ; in, and are the polynomial coefficients after long division. For any After z changes The rational form of is a discrete time sequence number, j is a positive integer index, is a positive integer power of n, for A step signal occurs at ; S43: Based on the sorted ,make , and obtain the shaping transfer function : ; in, is the pth numerator coefficient of the shaping transfer function, is the qth denominator coefficient of the shaping transfer function; S44: According to the linear system theory, using the numerator coefficient and denominator coefficient of the imaging transfer function, the recursive formula of the forming transfer function is obtained:

[0019] in, is the shaping transfer function The first denominator coefficient in is the shaping transfer function The N-1th denominator coefficient in is the shaping transfer function The first numerator coefficient in is the shaping transfer function The second numerator coefficient in is the shaping transfer function Middle The numerator coefficient, is the horizontal coordinate of the Nth control point.

[0020] Compared with the prior art, the invention can achieve the following beneficial effects: (1) The present invention creates the quasi-Gaussian pulse shaping method based on cubic spline interpolation, and simplifies the design and implementation process of the traditional shaping circuit system by using digital shaping design. The present invention is not limited by the constraints of fixed hardware circuit pulse shaping parameters and is not affected by the circuit environment. The present invention can achieve fast and stable pulse shaping by calculating the shaping transfer function through software such as MATLAB.

[0021] (2) The quasi-Gaussian pulse shaping method based on cubic spline interpolation created by the present invention has a smoother edge attenuation compared with the traditional trapezoidal shaping, which is due to the use of a fitted Gaussian shape; compared with other synthetic shaping methods, the cubic spline curve has the continuity of the second-order derivative, which can obtain a better low-pass filtering effect; compared with Salen-key shaping, the present invention has a clearer pulse width definition and more accurate pulse positioning and identification capabilities; compared with Gaussian window shaping, the pulse amplitude accuracy can be improved by adjusting the number of control points.

[0022] (3) The quasi-Gaussian pulse shaping method based on cubic spline interpolation created by the present invention can provide different accuracy solutions for different application scenarios. The present invention retains the advantages of the Gaussian shape, including the pulse sampling point width L and the Gaussian shape parameter α. It can appropriately increase the pulse shaping width and low-pass filtering capability in a low pulse count rate environment, thereby improving the pulse amplitude accuracy. By adjusting the Gaussian shape parameter α, the present invention can change the cutoff frequency and high-frequency attenuation coefficient of the photon counting imaging system in the frequency domain. By balancing the relationship between the cutoff frequency and the attenuation coefficient, the optimal parameter α is found to optimize the pulse amplitude accuracy.

[0023] (4) The quasi-Gaussian pulse shaping method based on cubic spline interpolation created by the present invention further improves the pulse amplitude accuracy by adding the shaping variable of the control point, including the number N of interpolation control points and the position of the interpolation control points. , … … The present invention can achieve effective approximation to Gaussian shaping by increasing the number of control points N. In addition, the present invention can also continuously try to adjust the positions of the control points while reducing the number of control points N, and fine-tune the shape while retaining the gentle rising edge and the gentle top of the Gaussian curve, and in the frequency domain of the transfer function of the photon counting imaging system, it can achieve a narrower low-pass transition band effect.

[0024] (5) The quasi-Gaussian pulse shaping method based on cubic spline interpolation created by the present invention is effectively adapted to the amplitude extraction of various exponential decay signals. The present invention fits the output curve of the pole-zero cancellation through a three-exponential model. Since the time decay constant of the original signal output by the photon counting imaging system fluctuates in the actual detection environment, the uncertainty is increased, resulting in over-compensation or under-compensation of the pole-zero cancellation. The time constant of the pole-zero cancellation signal can be further verified by fitting the three-exponential model. Compared with not using the exponential model fitting, better deconvolution impact shaping and better symmetric pulse shaping can be achieved.

[0025] (6) The quasi-Gaussian pulse shaping method based on cubic spline interpolation created by the present invention is universal and can be used to calculate the pulse amplitude of pulse sequences with exponential decay signals from different devices. The present invention can also be extended to non-Gaussian types and quadratic or higher-order interpolation fitting, such as Hamming window curves, Lorentz curves, and other customized bell-shaped curves with a middle bulge and attenuation on both sides. BRIEF DESCRIPTION OF THE DRAWINGS

[0026] The drawings constituting part of the present invention are used to provide a further understanding of the present invention. The exemplary embodiments and descriptions of the present invention are used to explain the present invention and do not constitute an improper limitation on the present invention. In the drawings: Figure 1 A schematic flow chart of a quasi-Gaussian pulse shaping method based on cubic spline interpolation according to an embodiment of the present invention; Figure 2 A schematic diagram of the structure of the photon counting imaging system according to an embodiment of the present invention; Figure 3 A schematic diagram of a signal output by a preamplifier according to an embodiment of the present invention; Figure 4 A schematic diagram of a single pulse signal output by a preamplifier according to an embodiment of the present invention; Figure 5 A schematic diagram of the change in the attenuation time constant of a single pulse signal output by a preamplifier according to an embodiment of the present invention; Figure 6 A schematic diagram of comparison between a preamplifier output signal and a pole-zero cancellation signal according to an embodiment of the present invention; Figure 7 A schematic diagram of fitting the three-exponential model described in the embodiment of the present invention; Figure 8 A schematic diagram of comparison between a preamplifier output signal and a deconvolution signal according to an embodiment of the present invention; Fig. 9 A schematic diagram of a standard Gaussian window according to an embodiment of the present invention; Fig.10 A two-dimensional graph with interpolation control points of [0, 18, 22, 26, 44] as described in the embodiment of the present invention; Fig.11 A two-dimensional graph with interpolation control points of [0, 3, 22, 41, 44] as described in the embodiment of the present invention; Fig.12 A schematic diagram of a cubic spline interpolation fitting curve according to an embodiment of the present invention; Fig.13 A schematic diagram showing the comparison between the fitting curves of the two control points described in the embodiment of the present invention and the standard Gaussian window; Fig.14 A schematic diagram of cubic spline interpolation fitting Gaussian shaping according to an embodiment of the present invention; Fig.15 A schematic diagram showing the comparison between the cubic spline interpolation fitting Gaussian shaping and the deconvolution signal described in the embodiment of the present invention; Fig.16 A schematic diagram comparing the cubic spline interpolation fitting Gaussian shaping described in the embodiment of the present invention and the standard Gaussian window.

[0027] Description of reference numerals: 1. Light source; 2. Detector; 3. Preamplifier; 4. Acquisition card; 5. Computer. DETAILED DESCRIPTION

[0028] In order to make the purpose, technical solution and advantages of the invention more clear, the invention is further described in detail below in conjunction with the accompanying drawings and specific embodiments. It should be understood that the specific embodiments described here are only used to explain the invention and do not constitute a limitation of the invention.

[0029] It should be noted that, in the absence of conflict, the embodiments of the present invention and the features in the embodiments may be combined with each other.

[0030] In the description of the present invention, it should be understood that the terms "center", "longitudinal", "lateral", "up", "down", "front", "back", "left", "right", "vertical", "horizontal", "top", "bottom", "inside", "outside" and the like indicate positions or positional relationships based on the positions or positional relationships shown in the accompanying drawings, and are only for the convenience of describing the present invention and simplifying the description, rather than indicating or implying that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation, and therefore cannot be understood as a limitation on the present invention. In addition, the terms "first", "second", etc. are only used for descriptive purposes, and cannot be understood as indicating or implying relative importance or implicitly indicating the number of technical features indicated. Thus, features defined as "first", "second", etc. may explicitly or implicitly include one or more of the features. In the description of the present invention, unless otherwise specified, "multiple" means two or more.

[0031] In the description of the invention, it should be noted that, unless otherwise clearly specified and limited, the terms "installation", "connection" and "connection" should be understood in a broad sense, for example, it can be a fixed connection, a detachable connection, or an integral connection; it can be a mechanical connection or an electrical connection; it can be a direct connection, or it can be indirectly connected through an intermediate medium, or it can be the internal communication of two components. For ordinary technicians in this field, the specific meanings of the above terms in the invention can be understood according to specific circumstances.

[0032] The present invention will be described in detail below with reference to the accompanying drawings and in combination with embodiments.

[0033] The present invention utilizes an acquisition card to perform analog-to-digital conversion on the output signal of a preamplifier, and then processes the data sequence, realizes quasi-Gaussian pulse shaping by combining pole-zero cancellation, three-exponential simulation fitting of pulse parameters and cubic spline interpolation fitting, and then extracts the pulse amplitude.

[0034] Specifically, the exponential decay signal output by the preamplifier is obtained by using an acquisition card, and the exponential decay signal is preprocessed by pole-zero cancellation to improve the pile-up and baseline unevenness caused by the undershoot and pulse width in the exponential decay signal. Then, the pole-zero cancellation signal fragment after pole-zero cancellation is fitted by a three-exponential model to obtain the corresponding time decay constant for deconvolution operation. Next, the pole-zero cancellation signal is solved based on the deconvolution system (deconvolution recursive formula) to obtain the deconvolution signal for pulse shaping and pulse identification and positioning. Next, the Gaussian shape is fitted by cubic spline interpolation, specifically: by selecting control points for the standard Gaussian window, and then performing cubic spline interpolation on the control points, the shaping transfer function of the cubic spline interpolation fitting Gaussian shape is obtained. Finally, the deconvolution signal is input into the shaping transfer function to obtain a quasi-Gaussian shaping signal. The signal is subsequently involved in the pulse amplitude extraction process.

[0035] like Figure 1 As shown, the present invention provides a quasi-Gaussian pulse shaping method based on cubic spline interpolation, which specifically includes the following steps: S1: Acquire an exponential decay signal, and perform pole-zero processing on the exponential decay signal to obtain a pole-zero cancellation signal; S11: Input signal It is represented as an exponential decay signal, and the Laplace transform of the exponential decay signal is performed: ; ; in, is an input signal in exponential decay form, is the time interval of the exponential decay signal from the pulse start time to the time when the pulse decays to 1 / e of the peak value, is a natural constant, t is time, s is a complex variable, is the exponential decay signal after Laplace transform; S12: Set the first time decay constant ,when The second time decay constant , using the first time decay constant and the second time decay constant to calculate the pole-zero cancellation transfer function in the s domain : ; Where s is a complex variable; S13: Convert the pole-zero cancellation transfer function to the time domain and discretize it to obtain the pole-zero cancellation recursive expression. Discretize to obtain an exponential decay signal , the exponential decay signal Substitute the pole-zero cancellation recursive expression to obtain the pole-zero cancellation signal : ; ; ; ; ; in, is the sampling period, , b, c and d are all intermediate variables for calculating the pole-zero cancellation signal and have no physical meaning. n is the serial number of the discrete time point.

[0036] Furthermore, the three-index model is: ; in, It is a three-exponential model, A is the first fitting parameter, B is the second fitting parameter, C is the third fitting parameter, D is the fourth fitting parameter, E is the fifth fitting parameter, F is the sixth fitting parameter, A, B, C ≥ 0, A+B+C=0, |A|>|B|>|C|.

[0037] S2: fitting the pole-zero cancellation fragment using the three-exponential model to obtain fitting parameters; recursively deducing the deconvolution system transfer function based on the fitting parameters to obtain the deconvolution recursive formula, and solving the pole-zero cancellation signal based on the deconvolution recursive formula to obtain the deconvolution signal; S21: Cancellation signal at pole zero Select the pole-zero cancellation fragment , the pole-zero cancellation fragments were fitted using a three-exponential model to obtain fitting parameters; The pole-zero cancellation segment is a complete, independent and non-pile-up pulse signal selected from the pole-zero cancellation signal, and the pulse peak amplitude of the pulse signal is at least 10 times of the noise standard deviation.

[0038] S22: Substitute the fitting parameters into the three-exponential model used in step S2 for Z transformation and organize them into Rational fractions of S23: recursively deducing the deconvolution system transfer function based on the rational fraction in step S22 to obtain a deconvolution recursive formula; The deconvolution system transfer function of the three-exponential model is: ; in, is a three-index model, is the impact function in the z domain; S24: Input the pole-zero cancellation signal into a deconvolution recursive formula for solving, and obtain a deconvolution signal.

[0039] The deconvolution recursive formula of the three-exponential model is:

[0040] ; ; ; ; ; ; ; ; ; ; in, , , , , , , , and are the coefficients of the deconvolution system transfer function expressing the three-exponential model and have no physical meaning. is the deconvolution signal.

[0041] S3: setting a standard Gaussian window, performing cubic spline interpolation fitting calculation on the standard Gaussian window, and designing an interpolation curve formula group of the Gaussian pulse shape; S31: Set the standard Gaussian window w: ; ; Among them, the length of the standard Gaussian window is L-1, L is the number of sampling points, L is an odd number, , is the width factor, is the standard deviation of the Gaussian curve, and n is the serial number of the discrete time point; S32: Select N sampling points from the L sampling points as interpolation control points, and the N sampling points are symmetrically selected with the sampling point n0=(L-1) / 2+1 as the center; S33: Fitting N interpolation control points based on the cubic sample interpolation method to obtain a fitting curve, and the N interpolation control points divide the fitting curve into N segments, and calculating the interpolation curve formula group of the N segments of the fitting curve: ; ;

[0042] ;

[0043] ; ; in, is the first fitting curve, is the second fitting curve, is the i-th segment fitting curve, is the N-1th segment fitting curve, is the Nth segment fitting curve, is the x coordinate of the ith control point among the N control points, The polynomial coefficients of the i-th curve are obtained by cubic spline fitting using control points.

[0044] S4: extrapolating the interpolation curve formula group to obtain a recursive formula of the forming transfer function; S41: discretize the interpolation curve formula group, and perform a neighboring term subtraction operation on the discretized interpolation curve formula group: ; ;

[0045] ;

[0046] ; ; in, is the first term of Gaussian pulse shaping fitted by cubic spline interpolation, The second term of Gaussian pulse shaping fitted by cubic spline interpolation, is the i-th term of Gaussian pulse shaping with cubic spline interpolation, is the N-1th term of Gaussian pulse shaping with cubic spline interpolation, is the Nth term of Gaussian pulse shaping with cubic spline interpolation, is the first segment of the cubic spline fitting curve after discretization, is the second segment of the cubic spline fitting curve after discretization, is the i-1th segment cubic spline fitting curve after discretization, is the discretized cubic spline fitting curve of the ith segment, is the N-2th segment of the discretized cubic spline fitting curve, is the N-1th segment of the discretized cubic spline fitting curve, For A step function where the position changes stepwise; S42: Use long division to To sort out and sort out Perform a Z transform: ; in, and are the polynomial coefficients after long division. For any After z changes The rational form of is a discrete time sequence number, j is a positive integer index, for A step signal occurs at ; S43: Based on the sorted ,make , and obtain the shaping transfer function : ; in, is the pth numerator coefficient of the shaping transfer function, is the qth denominator coefficient of the shaping transfer function; S44: According to the linear system theory, using the numerator coefficient and denominator coefficient of the imaging transfer function, the recursive formula of the forming transfer function is obtained:

[0047] in, is the shaping transfer function The first denominator coefficient in is the shaping transfer function The N-1th denominator coefficient in is the shaping transfer function The first numerator coefficient in is the shaping transfer function The second numerator coefficient in is the shaping transfer function Middle The numerator coefficient, is the horizontal coordinate of the Nth control point.

[0048] S5: Input the deconvolution signal into the shaping transfer function recursive formula and output a quasi-Gaussian shaping signal.

[0049] It should be noted that, in step S12, the first time decay constant of the pole-zero cancellation system is set ,when The pole-zero cancellation effect is achieved when the time constant As an estimate, the undershoot degree of the pole-zero cancellation output signal can be continuously adjusted. value until the undershoot phenomenon of the output pulse signal is effectively suppressed or eliminated. It is the artificially set pole-zero cancellation output signal attenuation time constant.

[0050] It should be noted that, in step S21, As a digital sequence, substitute n into the three-exponential model and make , and obtain the corresponding fitting parameters.

[0051] It should be noted that, in step S22, taking the three-index model as an example: The three-index model Performing Z transform yields: ; The formula Arrange into The rational fraction of : .

[0052] It should be noted that in step S32 and step S33, if , the N sampling points include n=-(L-1) / 2, n=0, and n=(L-1) / 2, which are the start and end points and the center point of the quasi-Gaussian shape.

[0053] The N sampling points are used as interpolation control points to perform cubic spline interpolation.

[0054] N control points p1, p2, ..., pN divide the entire fitting curve into N segments, that is to , and finally all 0 segment, among which, and The derivative of the outer endpoints p1 and pN of the curve segment is the endpoint gradient of the Gaussian curve.

[0055] It should be noted that in step S42, the known z-transformation formula is used to derive : ; ; ; ; ; in, Expressed as The rational form of .

[0056] It should be noted that the impulse function signal obtained in step 3 , as the input of the recursive formula of the forming transfer function, the expression is:

[0057] ;

[0058] ; The output obtained That is, cubic spline interpolation fits the Gaussian shaped signal.

[0059] Example 1 like Figure 2 As shown, a photon counting imaging system is constructed, and its process is as follows: light source 1 is used to make the incident light hit the microchannel plate, and after photomultiplication, a photoelectron group is output, and the photoelectron group falls on the two-dimensional position-sensitive anode of the detector 2, generating a pulse electrical signal with the position information of the photoelectron group, and the pulse electrical signal passes through the pre-charge sensitive amplifier (referred to as preamplifier 3) and is output to the subsequent shaping system to accurately identify and extract the peak value of the pulse signal. The photon counting imaging system subsequently uses the peak information of the pulse signal to accurately locate the position of the single photon, and then realizes photon counting imaging. The present invention only focuses on the relevant operations of shaping the digital pulses input to the computer 5 during the photon counting imaging process, and does not involve specific amplitude extraction methods and photon imaging calculations. In this imaging system, the acquisition card 4 can be used to sample the output analog signal of the preamplifier 3, and a cubic spline interpolation fitting Gaussian shaping experiment can be performed in the computer 5.

[0060] like Figure 3 As shown in the figure, the exponential decay signal output by the preamplifier 3 has obvious accumulation and undershoot, and the decay speed of each exponential signal is different, which will affect the baseline stability of the final formed signal. The DC signal in the pulse-free period is adjusted to 0 to achieve baseline zeroing. In the period with low pulse count rate, a representative non-pile-up pulse is selected for analysis. The representative non-pile-up pulse can be a single pulse signal with a small undershoot and a suitable decay speed, such as Figure 4 As shown, the decay rate of a single pulse signal is analyzed. The decay rate is the time length from the moment the signal occurs to the moment when the pulse height is 1 / e, as shown in Figure 5shown.

[0061] Step 1: Since the signal itself has noise and the decay time constant of each single pulse signal is different, only a rough estimate is needed, such as Figure 5 The estimated value of the pulse shown =130, then the pole-zero cancellation parameter of the pulse is =130, = 10. In order to eliminate the possibility of contingency, more representative non-pile-up pulses can be selected to estimate the decay time constant. According to experimental data, the decay time constant Pulses from 70 to 220 may exist. According to references and experience, Selecting a higher value can eliminate the undershoot of most pulses. Lower pulses will appear to overcompensate, but their effect on the baseline is smaller and The undercompensation caused by this is smaller. Therefore, a larger decay time constant can be selected, that is, =220, =10. The result of pole-zero cancellation is as follows Figure 6 As shown, the blue one is the pulse signal output by the preamplifier 3, and the red one is the pole-zero cancellation signal.

[0062] Step 2: After the pole-zero cancellation, the width of most single pulse signals is compressed, the accumulation is reduced, and the decay time constant of the pulse signal is reduced. At this time, you can select a suitable exponential signal and use the least squares method to perform exponential fitting on the target signal to obtain the decay time constant of the pulse signal after the pole-zero cancellation. =10, so theoretically the decay time constant is also close to 10, and the three-exponential model is used for fitting. The three-exponential signal to be fitted needs to include some undershoots and undershoot poles. The general length is set to arrive If the length is too long or too short, it is easy to cause overfitting or underfitting. An effective fit should make the difference between the fitting curve and the original single pulse signal less than or equal to the standard deviation of the noise.

[0063] The initial parameters of the fit are set to D = , E and F are set to be less than .like Figure 7 This is the result of fitting a pulse signal after pole-zero cancellation using a three-exponential model in Matlab. Goodness of fit R 2 =0.99 and above is considered to be an excellent fit. The coefficients of the deconvolution recursive formula are calculated using the decay time constant after the pole-zero cancellation obtained by fitting. The deconvolution signal is obtained through the deconvolution recursive formula, such as Figure 8As shown. The deconvolved impulse function signal can be used to locate the pulse in the subsequent amplitude extraction system, and the accumulation degree can be inferred by the pulse interval to achieve accumulation rejection. In the present invention, the deconvolved signal obtained by deconvolution is used as the input signal of the subsequent pulse Gaussian shaping system. The deconvolved signal input can ensure that the output signal of the shaping system has a fixed pulse width and a shape that conforms to the cubic spline interpolation fitting Gaussian shape designed in the present invention.

[0064] Step 3: Set a standard Gaussian window, perform cubic spline interpolation fitting calculation on the standard Gaussian window, and design an interpolation curve formula group for the Gaussian pulse shape.

[0065] Step 4: Calculate the interpolation curve formula group to obtain the forming transfer function.

[0066] The following is a detailed example. Under the constraints of L=45 and α=2.66, the standard Gaussian window is set as follows. Fig. 9 As shown in the figure, different control point positions are adjusted to compare the impulse response differences between the cubic spline interpolation Gaussian shaping method and the standard Gaussian window. The control points are guaranteed to cover the first and last control points and the center point of the standard Gaussian window. On this basis, control points are added symmetrically. For example, the control points [0, 18, 22, 26, 44] and the control points [0, 3, 22, 41, 44] are used. Figure 10-11 shown.

[0067] Perform cubic spline interpolation fitting on these five control points to obtain the interpolation fitting expressions of the five curves. Among them, the interpolation fitting curve of the control points [0, 3, 22, 41, 44] is as follows Fig.12 shown.

[0068] Compare the standard Gaussian window, [0, 18, 22, 26, 44] control point shaping and [0, 3, 22, 41, 44] control point shaping, the shaping effect is as follows Fig.13 In the time domain, the cubic spline fitting Gaussian shape retains the smooth and continuous time characteristics of the standard Gaussian, and has the same gain at the control point position. In addition, by adjusting the number and position of the control points, fine-tuning of different Gaussian shapes can be achieved.

[0069] Step 5: Input the impact function signal into the recursive formula of the forming system, and output the cubic spline forming, i.e., the quasi-Gaussian forming function. The result is as follows: Fig.14 shown. Fig.14 Some pulse signals in the image have undershoot and long tails, which is caused by the instability of the decay time constant of the detector 2 signal itself. The pulse position is located by deconvolution impact, and the judgment threshold of deconvolution positioning is set artificially. In this example, the threshold is set to Set to 0.03, that is, the deconvolution signal exceeding 0.03 is defined as the position of the pulse event. When the spacing between adjacent pulse signals is less than the shaping width L=45, pulse accumulation occurs, and the two accumulated pulses are excluded. After excluding the accumulated pulse events, the remaining non-pile-up pulses are peak extracted. The results show that the cubic spline interpolation fitting Gaussian shaping has higher peak accuracy than the standard Gaussian window. Fig.16 As shown in the figure, under the condition of forming width L=45, α=2.66, the combination of control points is traversed to obtain the peak error RMSE curve. Among them, the optimal error of cubic spline interpolation fitting Gaussian forming is less than 0.035, which has a lower peak error than the standard Gaussian window forming and better amplitude accuracy.

[0070] It should be understood that the various forms of processes shown above can be used to reorder, add or delete steps. For example, the steps described in the disclosure of the present invention can be performed in parallel, sequentially or in different orders, as long as the desired results of the technical solution disclosed in the present invention can be achieved, and this document does not limit this.

[0071] The above specific implementations do not constitute a limitation on the protection scope of the present invention. It should be understood by those skilled in the art that various modifications, combinations, sub-combinations and substitutions can be made according to design requirements and other factors. Any modification, equivalent substitution and improvement made within the spirit and principle of the present invention should be included in the protection scope of the present invention.

Claims

1. A quasi-Gaussian pulse shaping method based on cubic spline interpolation, characterized in that: The specific steps include: S1: Acquire an exponential decay signal, and perform pole-zero processing on the exponential decay signal to obtain a pole-zero cancellation signal; S2: Use the three-exponential model to fit the pole-zero cancellation fragment and obtain the fitting parameters; The transfer function of the deconvolution system is recursively deduced based on the fitting parameters to obtain a deconvolution recursive formula, and the pole-zero cancellation signal is solved based on the deconvolution recursive formula to obtain a deconvolution signal; S3: setting a standard Gaussian window, performing cubic spline interpolation fitting calculation on the standard Gaussian window, and designing an interpolation curve formula group of the Gaussian pulse shape; S4: extrapolating the interpolation curve formula group to obtain a recursive formula of the forming transfer function; S5: Input the deconvolution signal into the shaping transfer function recursive formula and output a quasi-Gaussian shaping signal.

2. The quasi-Gaussian pulse shaping method based on cubic spline interpolation according to claim 1, characterized in that: Step S1 specifically includes the following steps: S11: Input signal It is represented as an exponential decay signal, and the Laplace transform of the exponential decay signal is performed: ； ； in, is an input signal in exponential decay form, is the time interval of the exponential decay signal from the pulse start time to the time when the pulse decays to 1 / e of the peak value, is a natural constant, t is time, s is a complex variable, is the exponential decay signal after Laplace transform; S12: Set the first time decay constant , the second time decay constant , using the first time decay constant and the second time decay constant to calculate the pole-zero cancellation transfer function in the s domain : ； Where s is a complex variable; S13: Convert the pole-zero cancellation transfer function to the time domain and discretize it to obtain the pole-zero cancellation recursive expression. Discretize to obtain an exponential decay signal , the exponential decay signal Substitute the pole-zero cancellation recursive expression to obtain the pole-zero cancellation signal : ； ； ； ； ； in, is the sampling period, , b, c and d are all intermediate variables for calculating the pole-zero cancellation signal and have no physical meaning. n is the serial number of the discrete time point.

3. The quasi-Gaussian pulse shaping method based on cubic spline interpolation according to claim 1, characterized in that: The three-index model is: ； in, It is a three-exponential model, A is the first fitting parameter, B is the second fitting parameter, C is the third fitting parameter, D is the fourth fitting parameter, E is the fifth fitting parameter, F is the sixth fitting parameter, A, B, C ≥ 0, A+B+C=0, |A|>|B|>|C|.

4. The quasi-Gaussian pulse shaping method based on cubic spline interpolation according to claim 1, characterized in that: Step S2 specifically includes the following steps: S21: Cancellation signal at pole zero Select the pole-zero cancellation fragment , the pole-zero cancellation fragments were fitted using a three-exponential model to obtain fitting parameters; S22: Substitute the fitting parameters into the three-exponential model used in step S2 for Z transformation and organize them into Rational fractions of S23: recursively deducing the deconvolution system transfer function based on the rational fraction in step S22 to obtain a deconvolution recursive formula; S24: Input the pole-zero cancellation signal into a deconvolution recursive formula for solving, and obtain a deconvolution signal.

5. The quasi-Gaussian pulse shaping method based on cubic spline interpolation according to claim 4, characterized in that: The pole-zero cancellation segment is a complete, independent and non-pile-up pulse signal selected from the pole-zero cancellation signal, and the pulse peak amplitude of the pulse signal is at least 10 times of the noise standard deviation.

6. The quasi-Gaussian pulse shaping method based on cubic spline interpolation according to claim 4, characterized in that: The deconvolution system transfer function of the three-exponential model is: ； in, is a three-index model, is the impact function in the z domain; The deconvolution recursive formula of the three-exponential model is: ； ； ； ； ； ； ； ； ； ； in, , , , , , , , and are the coefficients of the deconvolution system transfer function expressing the three-exponential model and have no physical meaning. is the deconvolution signal.

7. The quasi-Gaussian pulse shaping method based on cubic spline interpolation according to claim 1, characterized in that: Step S3 specifically includes the following steps: S31: Set the standard Gaussian window w: ； ； Among them, the length of the standard Gaussian window is L-1, L is the number of sampling points, L is an odd number, , is the width factor, is the standard deviation of the Gaussian curve, and n is the serial number of the discrete time point; S32: among the L sampling points, N sampling points are selected as interpolation control points, and the N sampling points are selected symmetrically with the sampling point n0=(L-1) / 2+1 as the center; S33: Fitting N interpolation control points based on the cubic sample interpolation method to obtain a fitting curve, and the N interpolation control points divide the fitting curve into N segments, and calculating the interpolation curve formula group of the N segments of the fitting curve: ； ； ； ； ； in, is the first fitting curve, is the second fitting curve, is the i-th segment fitting curve, is the N-1th segment fitting curve, is the Nth segment fitting curve, is the x-coordinate of the ith control point among the N sampling points, The polynomial coefficients of the i-th curve are obtained by cubic spline fitting using control points.

8. The quasi-Gaussian pulse shaping method based on cubic spline interpolation according to claim 7, characterized in that: Step S4 specifically includes the following steps: S41: discretize the interpolation curve formula group, and perform a neighboring term subtraction operation on the discretized interpolation curve formula group: ； ； ； ； ； in, is the first term of Gaussian pulse shaping fitted by cubic spline interpolation, The second term of Gaussian pulse shaping fitted by cubic spline interpolation, is the i-th term of Gaussian pulse shaping with cubic spline interpolation, is the N-1th term of Gaussian pulse shaping with cubic spline interpolation, is the Nth term of Gaussian pulse shaping with cubic spline interpolation, is the first segment of the cubic spline fitting curve after discretization, is the second segment of the cubic spline fitting curve after discretization, is the i-1th segment cubic spline fitting curve after discretization, is the discretized cubic spline fitting curve of the ith segment, is the N-2th segment of the discretized cubic spline fitting curve, is the N-1th segment of the discretized cubic spline fitting curve, For A step function where the position changes stepwise; S42: Use long division to To sort out and sort out Perform a Z transform: ； in, and are the polynomial coefficients after long division. For any After z changes The rational form of is a discrete time sequence number, j is a positive integer index, is a positive integer power of n, for A step signal occurs at ; S43: Based on the sorted ,make , and obtain the shaping transfer function : ； in, is the pth numerator coefficient of the shaping transfer function, is the qth denominator coefficient of the shaping transfer function; S44: According to the linear system theory, using the numerator coefficient and denominator coefficient of the imaging transfer function, the recursive formula of the forming transfer function is obtained: ]。 in, is the shaping transfer function The first denominator coefficient in is the shaping transfer function The N-1th denominator coefficient in is the shaping transfer function The first numerator coefficient in is the shaping transfer function The second numerator coefficient in is the shaping transfer function Middle The numerator coefficient, is the horizontal coordinate of the Nth control point.

Images