Real-time correction method, correction effect analysis method and system based on EPID image

Abstract

The application discloses a real-time correction method and correction effect analysis method and system based on an EPID image, relates to the technical field of radiotherapy quality assurance, and comprises the following steps: automatically screening a dose image to obtain a cumulative dose image; positioning four corner points of a distorted field of an electronic portal imaging device (EPID) image to obtain a distorted quadrilateral; defining a target rectangle according to nominal field parameters input by a user; and realizing mapping of the distorted quadrilateral to the target rectangle by using a perspective transformation matrix to complete image correction. The application has high real-time performance, is compatible with FLASH radiotherapy, does not require any manual intervention in the whole process, can simultaneously correct in-plane rotation and X / Y axis tilt of the EPID, and makes up for the fact that existing EPID schemes have low automation degree, need manual intervention and cannot automatically correct perspective distortion, so that the analysis precision is seriously insufficient.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of radiotherapy quality assurance, in particular, a real-time correction method, correction effect analysis method and system based on EPID images. BACKGROUND

[0002] In radiotherapy (such as using a medical linear accelerator), the symmetry of the field is a key quality control (QA) indicator to ensure that the radiation dose is accurately delivered to the target area of the patient. The traditional QA method mainly relies on two technologies: film or water tank scanning. Among them:

[0003] Film QA process: the technician first places a special radiochromic film in the phantom and performs radiation exposure. After exposure, the film usually needs to go through scanning, scaling and other steps, and the whole process takes a long time.

[0004] Water tank QA process: this is the "gold standard" of dose measurement. This method uses a phantom filled with water (water tank), and a detector (such as an ionization chamber) is slowly moved in the water under the drive of a three-dimensional mechanical arm, and the dose distribution is scanned and measured point by point.

[0005] The above two traditional schemes have serious technical defects:

[0006] Time-consuming and cumbersome process: both film and water tank methods involve complex operations and long waiting or scanning times, and cannot provide real-time feedback.

[0007] Incompatible with FLASH radiotherapy: FLASH (ultra-high dose rate) radiotherapy is characterized by extremely short beam time (usually less than 1 second). The mechanical arm of the water tank cannot complete the scanning; the chemical reaction of the film also makes it impossible to provide real-time QA. Therefore, the traditional QA method is not suitable for FLASH radiotherapy.

[0008] Although there are technologies that attempt to use EPID (Electronic Portal Imaging Device, EPID, an electronic portal imaging device integrated into a linear accelerator for real-time imaging) for QA, there are:

[0009] Need for human intervention: when analyzing, the technician often needs to manually frame the ROI (Region of Interest, region of interest) and other field regions to be analyzed on the software.

[0010] Geometric correction model errors or omissions: Existing software can only (or requires manual correction) correct in-plane rotation caused by slight rotation of the gantry or EPD, limiting it to image distortion caused by EPD rotation within the plane—that is, a rectangular field of view becomes a "rotated rectangle," while its opposite sides remain parallel. More seriously, they are completely unable to automatically correct for X-axis or Y-axis tilt caused by the EPD detector plane not being perpendicular to the beam center axis. This tilt leads to perspective distortion, where a perfect rectangular field of view is projected as a trapezoid on the EPD. Analyzing symmetry based on this fundamentally introduces serious geometric errors.

[0011] Furthermore, existing EPD schemes lack validation and confidence feedback: existing systems are "black boxes." They analyze the acquired images "blindly," and even if the user's positioning is severely tilted, or a mechanical malfunction of the collimator causes an error in the field size, the system may still (based on error correction) arrive at a false positive conclusion of "pass." They cannot provide users with any quantitative indicators regarding the reliability of this analysis. Summary of the Invention

[0012] The purpose of this invention is to provide a real-time correction method, correction effect analysis method, and system based on EPID images. This addresses the problems of traditional QA methods being time-consuming, cumbersome, and completely unapplicable to FLASH radiotherapy scenarios; and the low automation level, manual intervention required, and inability to automatically correct fluoroscopic distortion in existing EPID analysis software, resulting in severely insufficient analysis accuracy. Furthermore, it solves the problem of existing EPID schemes lacking validation and confidence feedback, making it impossible to quantify their reliability.

[0013] The present invention solves the above problems through the following technical solution:

[0014] A real-time correction method based on EPD images, comprising:

[0015] Automatically filter dose images to obtain cumulative dose images;

[0016] Based on the cumulative dose image, the four corner points of the distorted field of the electronic field imaging device EPID are located to obtain the distorted quadrilateral;

[0017] The target rectangle is defined based on the nominal field parameters input by the user. The perspective transformation matrix is ​​used to map the distorted quadrilateral to the target rectangle, thus completing the image correction.

[0018] Furthermore, the automatic dose image screening includes: continuously acquiring the EPD image stream, automatically identifying the beam state using signal statistical features, and automatically identifying and extracting the cumulative dose image containing the radiation dose.

[0019] Furthermore, the continuous acquisition of the EPD image stream, the automatic identification of beam state using signal statistical features, and the automatic identification and extraction of the cumulative dose image containing the radiation dose include:

[0020] Maintain a first-in-first-out (FIFO) image frame buffer of length M in memory to cache the image stream;

[0021] N frames of images in the pre-beaming (waiting) state are acquired from the image frame buffer as the baseline, and the mean of the average pixel value (MPV) of the baseline image sequence is calculated. and standard deviation And the mean of the pixel standard deviation SD of the background image sequence. and standard deviation ;

[0022] Calculate the intensity threshold based on the statistical characteristics of the background. and noise threshold :

[0023] ;

[0024] ;

[0025] in, , This is the sensitivity coefficient;

[0026] Real-time calculation of the current frame The statistical characteristics of the data are such that when the formula (1) is satisfied, it is determined that the beam has started, and the data acquisition state is entered from the waiting state. The current frame and all subsequent frames are marked as dose image frames.

[0027] (1);

[0028] in, It is the average pixel value of the current frame. It is the standard deviation of pixels in the current frame;

[0029] In the data acquisition state, if continuous If the statistical characteristics of the frame image satisfy formula (2), then the beam is determined to be finished and acquisition is stopped;

[0030] (2);

[0031] Perform a pixel-by-pixel cumulative averaging operation on all dose image frames to generate a high signal-to-noise ratio cumulative dose image. .

[0032] Furthermore, the waiting state also includes dynamically updating the baseline, and the update method is as follows:

[0033] ;

[0034] in, This is the update rate coefficient; The mean of the average pixel values ​​(MPV) of the background image sequence before the update. The mean of the average pixel values ​​(MPV) of the updated background image sequence;

[0035] After dynamic baseline updates, adjust the intensity threshold. : .

[0036] Furthermore, the four corner points of the EPID image distortion field of the positioning electronic field imaging device include:

[0037] The penumbra region of interest (ROI) is extracted from the cumulative dose image based on the physical definition of isodose lines.

[0038] The distorted field boundary is located in the penumbra ROI, and straight line fitting and corner point calculation are performed on the distorted field boundary to obtain the field edge and four corner points.

[0039] Further, the extraction of the penumbra region of interest (ROI) from the cumulative dose image includes:

[0040] Traverse the cumulative dose image and search for the global maximum pixel value.

[0041] Based on the global maximum pixel value And the definition of the penumbra ROI, and calculation of the first isodose threshold. Second isodose threshold :

[0042] ;

[0043] ;

[0044] in, This represents the average value of the image edge region;

[0045] According to the first isodose threshold Second isodose threshold The cumulative dose image is binarized to obtain a binarized image. After removing isolated noise points from the binarized image, all closed contours are extracted, and the contours of the connected components with the largest area are retained to obtain the inner contour lines. and outer contour line ;

[0046] Inner contour line and outer contour line The enclosed ring-shaped region is the penumbra region (ROI).

[0047] Further, the step of locating the distorted field boundary of the penumbra ROI includes:

[0048] Step A1: Extract a 1D profile based on the penumbra ROI to obtain profile data, specifically including:

[0049] The penumbra ROI is sectioned along the horizontal / vertical direction, and each row / column covered by the penumbra is traversed to obtain 1D pixel value data. Based on this data, a 1D profile data point set is constructed. The profile data point set is a position-grayscale set. With position as the horizontal axis and grayscale as the vertical axis, an S-shaped curve is drawn.

[0050] Step A2: Use the nonlinear least squares method to fit the profile data to the fitting function. :

[0051] ;

[0052] in, It is the location index value in the profile data. These are all parameters that need to be obtained through fitting. The background value, For dosage range, The width of the penumbra, and That is, the x-coordinate of the center point of the S-curve;

[0053] Step A3: Calculate the values ​​for all rows / columns of the 1D profile using a fitting method. Parameters, to obtain The fitting results generate a set of sub-pixel point clouds for four edges. .

[0054] Furthermore, step A3 specifically includes:

[0055] Iterate through each row / column covered by the penumbra ROI;

[0056] Extract the pixel grayscale value sequence of the row / column within the penumbra ROI to construct a 1D profile data point set. ;in It is the index of the point's position in the grayscale value sequence. It is the gray value at the corresponding index in the gray value sequence; These are the corresponding row / column coordinates;

[0057] Define loss function : ,in It is the position sequence calculated based on the fitting function during the optimization process. The grayscale value result at the location;

[0058] Set the optimization goal as follows: convergence;

[0059] Setting parameters The initial value;

[0060] Perform nonlinear least squares fitting and iterative optimization. Four parameters, until The optimization stops when the difference from the actual pixel data is minimized, i.e., when the loss value is less than the set threshold or when the set maximum number of iterations is reached;

[0061] Converged parameters That is, the x-coordinate of the center point of the S-curve;

[0062] Repeat the above process for all rows / columns to generate a set of sub-pixel point clouds for the four edges. .

[0063] Furthermore, the process of performing line fitting and corner point calculation on the firing field boundary to obtain the firing field edge and four corner points includes:

[0064] Step B1: For the point cloud of each edge, perform the following steps respectively:

[0065] Step B11: Randomly select two points from the point cloud and calculate a straight line model;

[0066] Step B12: Traverse all points in the point cloud and calculate the algebraic distance from each point to a line. If the distance is less than a preset threshold... Mark them as interior points;

[0067] Step B13: Repeat steps B11 to B12 a set number of times, and retain the line model with the most interior points as the best coarse model.

[0068] Step B14: Extract all interior points corresponding to the best coarse model, and use the least squares method to perform linear regression to obtain the final high-precision straight line equation;

[0069] Step B2: Solve the four fitted high-precision line equations simultaneously to obtain the sub-pixel coordinates of the four intersection points. Simultaneously calculate the lengths of all four sides.

[0070] Further, the step of defining a target rectangle based on the nominal field parameters input by the user, and using a perspective transformation matrix to map the distorted quadrilateral to the target rectangle to complete image correction includes:

[0071] Receive user input of nominal field physical width and nominal field physical height And convert it to pixel units:

[0072] ;

[0073] ;

[0074] in, The pixel-to-millimeter conversion factor for EPD; The nominal pixel width of the field of view; The nominal pixel height of the shooting field;

[0075] Determine the coordinates of the four vertices of the target rectangle, and set the height and width of the target rectangle to be equal to... and ;

[0076] Using perspective transformation matrix This associates the coordinates of points on the distorted image with the corresponding coordinates of points on the target rectangle.

[0077] Solve for the perspective transformation matrix using the four vertices of the target rectangle. This method achieves a precise mapping from a distorted quadrilateral to a target rectangle. Specifically, it substitutes the coordinates of four corresponding points to construct a system of linear equations, and then solves the system of equations using singular value decomposition (SVD) to obtain the normalized perspective transformation matrix. This achieves a precise mapping from the distorted quadrilateral to the target rectangle, resulting in the target image.

[0078] Traverse the coordinates of each pixel in the target image Using matrices inverse matrix Mapping the target image coordinates back to the source image, i.e., the distorted quadrilateral coordinate system, yields floating-point coordinates. ;

[0079] Select floating-point coordinates in the distorted quadrilateral Surrounding Using a neighborhood of pixels, a cubic convolution kernel function is applied. Calculate the weighted average gray value and assign it to the pixel coordinates. ;

[0080] Generate a geometrically corrected dose image .

[0081] Furthermore, the perspective transformation matrix for:

[0082] ;

[0083] in, For affine components, used to correct image rotation caused by collimator rotation or rotation in the EPD plane, as well as overall scaling caused by source-image distance deviation; The translation component is used to correct the X-axis and Y-axis translation deviations of the field center relative to the imaging center; The perspective component is used to correct the trapezoidal distortion caused by the EPD detector plane not being perpendicular to the beam center axis; These are the normalized parameters.

[0084] A method for analyzing the correction effect based on EPD images includes obtaining a corrected image using the real-time correction method based on EPD images, calculating symmetry based on the corrected image, and calculating confidence based on the difference between the distorted quadrilateral and the corrected image.

[0085] Furthermore, the calculation of symmetry based on the corrected image includes:

[0086] On the corrected image, the horizontal dose distribution curve and the vertical dose distribution curve are extracted along the horizontal and vertical central axes, respectively.

[0087] Calculate the horizontal and vertical symmetry within the field width n3:

[0088] Horizontal symmetry is the maximum of the absolute values ​​of the ratio of the difference between symmetrical pixels on the horizontal dose distribution curve to the maximum value on the horizontal dose distribution curve.

[0089] Vertical symmetry is the maximum of the absolute values ​​of the ratio of the difference between symmetrical pixels on the vertical dose distribution curve to the maximum value on the vertical dose distribution curve.

[0090] Furthermore, the confidence score calculation based on the difference between the distorted quadrilateral and the corrected image includes:

[0091] Calculate the confidence level for geometric distortion, including calculating the Y-axis tilt distortion factor and the X-axis tilt distortion factor for trapezoidal distortion, where:

[0092] The Y-axis tilt distortion factor is the absolute value of the difference between the ratio of the upper side length to the lower side length of the distorted quadrilateral and 1.

[0093] The X-axis tilt distortion factor is the absolute value of the difference between the ratio of the left length of the distorted quadrilateral to the right length of the distorted quadrilateral and 1.

[0094] The confidence level of dimensional fidelity is calculated based on comparing the measured size with the nominal size, including:

[0095] Calculate the average width and average height of the distorted quadrilateral:

[0096] The average width is half the sum of the upper and lower sides of the distorted quadrilateral.

[0097] The average height is half the sum of the left and right lengths of the distorted quadrilateral.

[0098] Calculate the percentage of width and height scale errors between the distorted quadrilateral and the corrected image:

[0099] The width scale error percentage is the absolute value of the ratio of the difference between the average width and the nominal pixel width to the nominal pixel width;

[0100] The height scale error percentage is the absolute value of the ratio of the difference between the average height and the nominal pixel height to the nominal pixel height;

[0101] The nominal pixel width is the product of the nominal physical width of the field of view input by the user and the pixel-to-millimeter conversion factor of EPD, and the nominal pixel height is the product of the nominal physical length of the field of view input by the user and the pixel-to-millimeter conversion factor of EPD.

[0102] Furthermore, it also includes visualizing the symmetry calculation results and confidence calculation results, and generating reports.

[0103] A real-time automatic correction system for radiotherapy radiation fields based on EPID images and perspective transformation, comprising:

[0104] The user input module is configured to allow users to input nominal field parameters;

[0105] The dose image filtering module is configured to automatically filter dose images to obtain cumulative dose images;

[0106] The field distortion parameter calculation module is configured to locate the four corner points of the distorted field of the electronic field imaging device EPID image based on the cumulative dose image, and obtain the distorted quadrilateral.

[0107] The perspective transformation correction module is configured to define the target rectangle based on the nominal field parameters input by the user, and use the perspective transformation matrix to map the distorted quadrilateral to the target rectangle, thereby completing the image correction.

[0108] Furthermore, it also includes:

[0109] The analysis and confidence calculation module is configured to calculate symmetry based on the corrected image and confidence based on the difference between the distorted quadrilateral and the corrected image;

[0110] The visualization and reporting module is configured to visualize the symmetry calculation results and confidence calculation results and generate reports.

[0111] Compared with the prior art, the present invention has the following advantages and beneficial effects:

[0112] (1) This invention utilizes the real-time imaging characteristics of EPD to complete the analysis within seconds after exposure. It has high real-time performance and is perfectly compatible with the ultra-short beam exit time of FLASH radiotherapy. It solves the problems of long time consumption, cumbersome operation and inability to be applied to FLASH radiotherapy scenarios of traditional water tank and film methods.

[0113] (2) The present invention automatically completes the shooting field positioning, distortion detection and perspective correction without any manual intervention, eliminating the subjective error caused by manual selection of ROI and other operations.

[0114] (3) This invention employs perspective transformation, which can simultaneously correct in-plane rotation and X / Y axis tilt (perspective distortion) of EPD, thus overcoming the shortcomings of existing EPD schemes, such as low automation, the need for manual intervention, and the inability to automatically correct perspective distortion, which leads to severely insufficient analysis accuracy. More importantly, this invention uses user-input nominal parameters (rather than estimated values) as the correction benchmark, fundamentally eliminating scale estimation errors in the correction process and ensuring the geometric accuracy and scale accuracy of the corrected image.

[0115] (4) By introducing the nominal field parameters input by the user as a benchmark, the present invention realizes the multi-dimensional confidence assessment of the analysis results, which solves the defects of blind analysis in the prior art; it provides decision assistance to users and prevents QA misjudgment caused by user input errors or equipment mechanical failures.

[0116] (5) This invention can quantitatively assess the reliability of the EPD's positioning posture through geometric distortion confidence. It can also quantitatively assess the consistency between the measurement scale and the nominal scale through scale fidelity confidence. Attached Figure Description

[0117] Figure 1 This is a flowchart of an embodiment of the present invention;

[0118] Figure 2 This is a system principle block diagram according to an embodiment of the present invention;

[0119] Figure 3This is a schematic diagram of the calculation of field distortion parameters in an embodiment of the present invention, wherein (a) is a schematic diagram of the extraction of the penumbra ROI; (b) is a schematic diagram of the fitting curve; and (c) is a schematic diagram of the edge sub-pixel point cloud.

[0120] Figure 4 This is a schematic diagram of perspective correction based on the nominal field of view according to an embodiment of the present invention;

[0121] Figure 5 This is a schematic diagram of the system display interface according to an embodiment of the present invention. Detailed Implementation

[0122] The present invention will be further described in detail below with reference to embodiments, but the implementation of the present invention is not limited thereto.

[0123] Example 1:

[0124] Combined with appendix Figure 1 As shown, a real-time correction method based on EPD images includes:

[0125] (a) User inputs nominal field parameters

[0126] The input interface allows users to input the nominal field parameters of QA, including the nominal field physical width. and nominal field physical height (Units are usually millimeters (mm); the pixel-to-millimeter conversion factor of the EPD needs to be pre-calibrated. The unit is pixels per millimeter. Alternatively, the coefficient can be automatically calculated when the source-detector distance (SID) changes for later use.

[0127] (ii) Automatically filter dose images to obtain cumulative dose images

[0128] The image stream continuously acquired from the EPD automatically identifies the beam state using signal statistical characteristics, and automatically identifies and extracts dose images containing radiation dose. This includes:

[0129] Due to the complex environment of the accelerator, random electromagnetic interference pulses may exist. If triggered by a single frame image, misjudgment is likely to occur. Therefore, in this embodiment, the image stream continuously acquired from the EPD is first buffered in memory through a sliding window into an image frame buffer (first-in-first-out FIFO queue with a length of M, such as M=3). Only when M consecutive frames of data show stable dose characteristics is it confirmed as a real beam output.

[0130] Automatic acquisition of beam before Frames (e.g.) Using the background image as the baseline, background statistical modeling is achieved. The mean of the average pixel value (MPV) of the background image sequence is calculated. and standard deviation and the mean of pixel standard deviation (SD) and standard deviation .

[0131] Based on the background (the background being an image with no dose, where grayscale values ​​are uniformly distributed throughout the image) statistical characteristics, the trigger threshold is automatically calculated:

[0132] Intensity threshold : ;

[0133] noise threshold :

[0134] in, This is the sensitivity coefficient, for example, 5.

[0135] Real-time calculation of the current frame Statistical characteristics:

[0136] 1) When the buffer contains consecutive Frame image does not satisfy: and At this time, the system is in a "waiting" state and continues to update. This is to accommodate zero-point drift caused by temperature changes or prolonged operation of the detector. Zero-point drift is mainly manifested as an overall shift in the reference voltage of all pixels, i.e., the entire image... Figure 1 The brightness or darkness changes slightly, so this drift primarily alters the MPV (MPV) and has minimal impact on the SD (SD). The update formula is:

[0137] ;

[0138] in, The update rate coefficient (typically takes the value of ) ).renew Then, the intensity threshold is adjusted according to the formula below: This enables dynamic baseline updates.

[0139] 2) When the buffer contains consecutive All frames of images satisfy:

[0140] ;

[0141] Once the beam emission is detected, the system status changes from "waiting" to "acquisition," and the current frame and all subsequent frames are marked as dose images; among them, It is the average pixel value of the current frame. It is the standard deviation of pixels in the current frame.

[0142] 3) In the data acquisition state, when continuously... Frames (e.g.) The statistical characteristics of the image satisfy:

[0143] or If the signal is lost, the beam is considered to have ended, and data collection is stopped.

[0144] After beam output, a pixel-by-pixel cumulative averaging operation is performed on all selected dose image frames to generate a high signal-to-noise ratio cumulative dose image. .

[0145] (iii) Based on the cumulative dose image, locate the four corner points of the distorted field of the electronic field imaging device EPID image to obtain the distorted quadrilateral.

[0146] To overcome the interference from EPD image noise and penumbra, a two-stage, coarse-to-fine localization method is adopted, and the four corner points of the firing field are calculated with sub-pixel accuracy. This embodiment does not directly use a simple pixel threshold to determine edges, but instead employs a "mathematical model fitting" method to achieve sub-pixel-level accurate measurement of the firing field edges. For example... Figure 3 As shown, it specifically includes:

[0147] 1. Extraction of the penumbra region ROI based on the physical definition of isodose lines

[0148] (1) Dose normalization: Traversing the cumulative dose image Search for the global maximum pixel value (Representing 100% dose point) and the average value of the image edge region (Represents 0% background);

[0149] (2) Calculation of physical threshold: Calculate the threshold of the first isodose line according to the radiation physics industry standard (definition of 20%-80% penumbra). Second isodose threshold :

[0150] ;

[0151] .

[0152] (3) Topological contour extraction: Apply a threshold to the image and Binarization is performed on each image. Morphological opening (erosion followed by dilation) is then performed on the binary images to remove isolated noise. A contour tracking algorithm (such as the Suzuki85 algorithm) is used to extract all closed contours, retaining the contours of the connected components with the largest area, thus obtaining the inner contour lines. and outer contour line .

[0153] (4) Defining the penumbra ROI: The penumbra ROI is strictly defined as and The ring-shaped area enclosed by these two contour lines (such as...) Figure 3 (As shown in (a)). Subsequent high-precision calculations are performed only within this region, shielding the high-dose flat area and the background noise area.

[0154] 2. Sub-pixel precise positioning of the firing field boundary

[0155] Sub-pixel accuracy refers to the use of mathematical methods in digital image processing to improve the precision of measurement or positioning to smaller units (such as 1 / 2 or 1 / 4 of a pixel) between physical pixels, thereby breaking through the limitations of hardware resolution and achieving higher precision calculation and analysis.

[0156] The physical pixels of EPID are discrete (e.g., one pixel represents 0.4mm). If the physical edge of the field of view falls exactly between the 100th and 101st pixels (e.g., at position 100.4), traditional pixel-level methods can only determine whether it is 100 or 101, resulting in an error of 0.4mm. This embodiment achieves sub-pixel accuracy through function fitting, which can calculate the floating-point coordinate "100.2".

[0157] (1) Extracting 1D profile

[0158] For the left / right edges of the penumbra ROI, section along the horizontal direction (extract along the X-axis); for the top / bottom edges, section along the vertical direction (extract along the Y-axis), and traverse each row (for the left and right edges) or each column (for the top and bottom edges) covered by the penumbra. The traversal results in 1D pixel value data, and a profile data point set is constructed based on this data.

[0159] left edge For example, iterate through each row covered by the "left penumbra ROI". ,here These are the row coordinates of the image. Extract the pixel grayscale value sequence within the penumbra ROI of that row to construct a 1D profile data point set. ,in It is the index of the point's position in the grayscale value sequence. It is the gray value at the corresponding index in the gray value sequence.

[0160] (2) Nonlinear least squares fitting

[0161] Because the penumbra physically exhibits an S-shaped (Sigmoid) gradient, this invention uses a nonlinear least squares method to fit the 1D profile data (which is essentially a "position-grayscale" set; if position is plotted on the horizontal axis and grayscale on the vertical axis, the curve resembles the letter S, closely approximating the shape of the Sigmoid function) to a Sigmoid function model. This invention constructs the following Sigmoid fitting function:

[0162] ;

[0163] Where x is the location index value in the profile data. These are all parameters that need to be obtained through fitting. The background value, For dosage range, The width of the penumbra, and This is the x-coordinate corresponding to the center point of the S-curve. It is calculated using a fitting method. The parameters are used to enable the algorithm to find a curve that best fits the current situation based on the specific pixel data of each row, thereby obtaining a more accurate result. The fitting results, such as Figure 3 As shown in (b).

[0164] Here, they must pass through simultaneously, either in rows or columns. and These two contour lines are isodose lines.

[0165] When the fit is optimal, the parameters Mathematically, this corresponds precisely to the 50% center position of the dose change. The mathematical principle is that when a point is taken on the horizontal axis... , making Exactly equal to the parameter Time (i.e.) The following results can be obtained:

[0166] ;

[0167] in, This is precisely "the background value plus half the amplitude," which is the 50% dose point, and the 50% dose point is the definition of the edge of the radiation field. The detailed fitting process is as follows:

[0168] a. Define the loss function (sum of squared residuals): .in It is the position sequence calculated based on the fitting function during the optimization process. The output is grayscale value. The optimization goal is to adjust... Parameters that make and When the values ​​are close, it means the Loss converges.

[0169] b. Setting initial values ​​for optimization parameters: To ensure fast convergence and high accuracy of the fit, the initial values ​​of the four parameters are first estimated based on the original data sequence. The minimum gray value in the profile data is set to... The initial value; set the dynamic range (maximum value minus minimum value) of grayscale values ​​in the profile data to . The initial value; the sequence values ​​in the profile data are closest to the median value ( The position of ) is set to The initial value; the pixel width of the upper penumbra ROI in the row (or column) where the profile data is located is used as the initial value. The initial value.

[0170] c. The Levenberg-Marquardt (LM) algorithm is used for nonlinear least-squares fitting. The algorithm iteratively optimizes... Four parameters, until the curve The optimization process stops when the difference from the actual pixel data is minimized, that is, when the loss value is less than a set threshold (typically 0.01). Furthermore, the optimization process also stops when the set maximum number of iterations is reached to avoid an infinite loop.

[0171] d. Parameters after convergence That is, the x-coordinate of the center point of the S-curve.

[0172] (3) Subpixel cloud

[0173] Repeat the above process for all rows / columns to generate a set of sub-pixel point clouds for the four edges. ,like Figure 3 As shown in (c).

[0174] 3. RANSAC line fitting and corner point calculation

[0175] Random Sample Consensus (RANSAC) is an algorithm that calculates mathematical model parameters of a dataset containing outliers to obtain valid sample data.

[0176] (1) For the point cloud of each edge (e.g.) Perform the following steps respectively:

[0177] a. Random sampling: Randomly select two points from the point cloud and calculate a straight line model.

[0178] b. Interior Point Statistics: Traverse all points in the point cloud and calculate the algebraic distance from each point to a line. If the distance is less than a preset threshold... (e.g., 0.5 pixels), marked as "Inlier".

[0179] c. Model Evaluation: Repeat the above steps. Next, the straight line model with the most "interior points" is retained as the best coarse model.

[0180] d. Final Refinement: Extract all "interior points" corresponding to the best model, perform linear regression using the least squares method, and obtain the final high-precision linear equation (e.g., ).

[0181] (2) Corner point solution: Solve the four fitted lines simultaneously ( The system of equations yields the sub-pixel coordinates of the four intersection points. Calculate the lengths of all four sides simultaneously. The length of the upper side of the distorted quadrilateral The length of the lower side of the distorted quadrilateral The left side of the distorted quadrilateral is longer. The right side of the distorted quadrilateral is longer. .

[0182] (iv) Define the target rectangle according to the nominal field parameters input by the user, and use the perspective transformation matrix to realize the mapping from the distorted quadrilateral to the target rectangle to complete the image correction.

[0183] This part implements high-precision correction of distorted field images into standard rectangles based on nominal parameters input by the user.

[0184] (1) Define the target space, such as Figure 4 As shown, it includes:

[0185] a. Convert nominal parameters to pixel units:

[0186] ; .

[0187] in, The nominal physical width of the firing field. The nominal physical height of the firing field; The pixel-to-millimeter conversion factor for EPD; The nominal pixel width of the field of view; The nominal pixel height of the shooting field.

[0188] b. Define the "ideal target rectangle" The coordinates of the four vertices Make its center located at the center of the image, with length and width respectively. and .

[0189] (2) Constructing the Homography Matrix

[0190] a. Define and solve a perspective transformation matrix This matrix defines the projection mapping relationship from the distorted image plane to the ideal standard plane.

[0191] Assume the coordinates of a point on the source image (distorted quadrilateral) are... The coordinates of the corresponding point on the target image (ideal target rectangle) are: Both are determined through the homography matrix. The correlation is shown in the following formula:

[0192] ;

[0193] Each parameter in the homography matrix has a clear geometric correction physical meaning:

[0194] (Affine component): Primarily corrects image rotation caused by collimator rotation or rotation in the EPD plane, as well as overall scaling caused by source-image distance (SID) deviation.

[0195] (Translation component): Corrects the X-axis and Y-axis translation deviations of the field center relative to the imaging center.

[0196] (Perspective components): They are responsible for correcting trapezoidal distortion caused by the EPD detector plane not being perpendicular to the beam center axis (i.e., tilting).

[0197] : Normalization parameter, usually constrained as follows .

[0198] b. Substitute the coordinates of the four corresponding points to construct a structure of the form shown below. The system of linear equations is obtained by solving the system of equations using singular value decomposition (SVD) to obtain the normalized equations. Perspective transformation matrix This matrix can realize the transformation from the source quadrilateral. (Depend on (Define) to target rectangle (Depend on The precise mapping of (definition).

[0199] (3) Inverse mapping and grayscale interpolation

[0200] a. Reverse scanning: Traversing the coordinates of every pixel in the target image (corrected image). Using matrices inverse matrix Mapping the target coordinates back to the source image coordinate system yields floating-point coordinates. .

[0201] b. Bicubic Interpolation: Selects from the source image... Surrounding Using a neighborhood of pixels, a cubic convolution kernel function is applied. Calculate the weighted average gray value and assign it to the target pixel. This interpolation method can preserve the smoothness and gradient information of the dose distribution to the greatest extent.

[0202] Generate a geometrically corrected dose image And output it.

[0203] Example 2:

[0204] This invention provides a method for analyzing the correction effect based on EPD images. The method utilizes the aforementioned real-time correction method based on EPD images to obtain a corrected image, calculates the symmetry based on the corrected image, and calculates the confidence level based on the difference between the distorted quadrilateral and the corrected image. The method includes:

[0205] (I) Symmetry Analysis

[0206] 1. Extract the horizontal dose distribution curves along the horizontal and vertical central axes from the corrected image, respectively. and vertical dose distribution curve ;

[0207] 2. Based on the IEC standard formula, calculate the symmetry within the width range of the firing field, preferably, for example, within 80% of the width range. :

[0208] Horizontal symmetry is the maximum of the absolute values ​​of the ratio of the difference between symmetrical pixels on the horizontal dose distribution curve to the maximum value on the horizontal dose distribution curve.

[0209] Vertical symmetry is the maximum of the absolute values ​​of the ratio of the difference between symmetrical pixels on the vertical dose distribution curve to the maximum value on the vertical dose distribution curve.

[0210] Taking horizontal symmetry as an example, the calculation formula is:

[0211] ;

[0212] in, The maximum value on the horizontal dose distribution curve. These are symmetrical pixels.

[0213] (II) Calculation of multi-dimensional confidence

[0214] 1. Confidence level of geometric distortion It is used to evaluate the placement.

[0215] Geometric distortion confidence Includes: Y-axis tilt distortion factor for calculating trapezoidal distortion and X-axis tilt distortion factor , Used to calculate the difference between the upper and lower bases. Used to calculate the difference between the left and right sides:

[0216] ;

[0217] ;

[0218] in, Let be the length of the top side of the distorted quadrilateral; Let be the length of the lower side of the distorted quadrilateral; The left side of the distorted quadrilateral is longer; The right side length of the distorted quadrilateral;

[0219] 2. Scale fidelity confidence level is used to evaluate the accuracy of the input field size / mechanical settings (such as collimator parameter settings). This is because when the confidence level value is too high, there are two possibilities: one is that the nominal field size input by the user does not match the actual size; the other is that the collimator parameters are set incorrectly, resulting in a field size that does not match the actual size.

[0220] This embodiment provides the following evaluation method:

[0221] (1) Compare the measured dimensions with the nominal dimensions

[0222] Calculate the average width of the distorted quadrilateral (pixels) and average height (pixels):

[0223] ;

[0224] ;

[0225] Get nominal width / height (pixels): and ;

[0226] (2) Calculate the percentage error in width scale between the distorted quadrilateral and the corrected image. and height scale error percentage :

[0227] ;

[0228] .

[0229] Furthermore, it also includes presenting the analysis results and multi-dimensional confidence indicators to users in an intuitive and easy-to-understand manner, and providing intelligent decision-making assistance. The displayed content includes:

[0230] The corrected image, such as Figure 5 Section A, the main display area, shows the calibrated EPID image and overlays X / Y analysis lines;

[0231] Contour curves: Plot one-dimensional dose distribution curves along the X and Y axes, such as... Figure 5 Section B. Contour curve area;

[0232] Symmetrical results, such as Figure 5 Section C. Symmetry Analysis Results.

[0233] Multi-dimensional confidence, such as Figure 5 Section D. Multidimensional Confidence Report (Decision Support) includes:

[0234] a. Geometric distortion: Display and The value is set to [value]. If it exceeds the threshold (e.g., >2%), the message "Warning: Severe trapezoidal distortion detected! Please check EPD placement" will be displayed.

[0235] b. Scale fidelity: Display and The value. If it exceeds the threshold (e.g., >3%), the following message will be displayed: "Warning: The detected field size is significantly different from the input value! Please check the user input or collimator settings." Figure 5 As shown.

[0236] Example 3:

[0237] Combination Figure 2 As shown, a real-time automatic correction system for radiotherapy radiation fields based on EPD images and perspective transformation includes:

[0238] The user input module is configured to allow users to input nominal field parameters;

[0239] The dose image filtering module is configured to automatically filter dose images to obtain cumulative dose images;

[0240] The field distortion parameter calculation module is configured to locate the four corner points of the distorted field of the EPID image of the electronic field imaging device, and obtain the distorted quadrilateral.

[0241] The perspective transformation correction module is configured to define the target rectangle based on the nominal field parameters input by the user, and use the perspective transformation matrix to map the distorted quadrilateral to the target rectangle, thereby completing the image correction.

[0242] Furthermore, it also includes:

[0243] The analysis and confidence calculation module is configured to calculate symmetry based on the corrected image and confidence based on the difference between the distorted quadrilateral and the corrected image;

[0244] The visualization and reporting module is configured to visualize the symmetry calculation results and confidence calculation results and generate reports.

[0245] Although the present invention has been described herein with reference to illustrative embodiments, the above embodiments are merely preferred embodiments of the present invention, and the implementation of the present invention is not limited to the above embodiments. It should be understood that those skilled in the art can devise many other modifications and implementations, which will fall within the scope and spirit of the principles disclosed in this application.

Claims

1. A real-time correction method based on EPD images, characterized in that, include: Automatically filter dose images to obtain cumulative dose images; Based on the cumulative dose image, the four corner points of the distorted field of the electronic field imaging device EPID are located to obtain the distorted quadrilateral; The four corner points of the EPID image distortion field of the positioning electronic field imaging device include: Extract the penumbra region of interest (ROI) from the cumulative dose image; Locate the distorted radiation field boundary in the penumbra region ROI; By performing line fitting and corner point calculation on the distorted field boundary, the field edge and four corner points are obtained. The process of locating the distorted field boundary of the penumbra ROI includes: Step A1: Extract a 1D profile based on the penumbra ROI to obtain profile data, specifically including: The penumbra ROI is sectioned along the horizontal / vertical direction, and each row / column covered by the penumbra is traversed to obtain 1D pixel value data. Based on this data, a 1D profile data point set is constructed. The profile data point set is a position-grayscale set. With position as the horizontal axis and grayscale as the vertical axis, an S-shaped curve is drawn. Step A2: Use the nonlinear least squares method to fit the profile data to the fitting function. : ； in, It is the location index value in the profile data. These are all parameters that need to be obtained through fitting. The background value, For dosage range, The width of the penumbra, and That is, the x-coordinate of the center point of the S-curve; Step A3: Calculate the values ​​for all rows / columns of the 1D profile using a fitting method. Parameters, to obtain The fitting results generate a set of sub-pixel point clouds for four edges. ; The process of fitting straight lines and calculating corner points on the distorted field boundary yields the field edge and four corner points, including: Step B1: For the point cloud of each edge, perform the following steps respectively: Step B11: Randomly select two points from the point cloud and calculate a straight line model; Step B12: Traverse all points in the point cloud and calculate the algebraic distance from each point to a line. If the distance is less than a preset threshold... Mark them as interior points; Step B13: Repeat steps B11 to B12 a set number of times, and retain the line model with the most interior points as the best coarse model. Step B14: Extract all interior points corresponding to the best coarse model, and use the least squares method to perform linear regression to obtain the final high-precision straight line equation; Step B2: Solve the four fitted high-precision line equations simultaneously to obtain the sub-pixel coordinates of the four intersection points. Simultaneously calculate the lengths of all four sides; The target rectangle is defined based on the nominal field parameters input by the user. The perspective transformation matrix is ​​used to map the distorted quadrilateral to the target rectangle, thus completing the image correction.

2. The real-time correction method based on EPD images according to claim 1, characterized in that, Extracting the penumbra region of interest (ROI) from the cumulative dose image includes: Traverse the cumulative dose image and search for the global maximum pixel value. ; Based on the global maximum pixel value And the definition of the penumbra ROI, and calculation of the first isodose threshold. Second isodose threshold : ； ； in, This represents the average value of the image edge region; According to the first isodose threshold Second isodose threshold The cumulative dose image is binarized to obtain a binarized image. After removing isolated noise points from the binarized image, all closed contours are extracted, and the contours of the connected components with the largest area are retained to obtain the inner contour lines. and outer contour line ; Inner contour line and outer contour line The enclosed ring-shaped region is the penumbra region (ROI).

3. The real-time correction method based on EPD images according to claim 1, characterized in that, Step A3 specifically includes: Iterate through each row / column covered by the penumbra ROI; Extract the pixel grayscale value sequence of the row / column within the penumbra ROI to construct a 1D profile data point set. ;in It is the index of the point's position in the grayscale value sequence. It is the gray value at the corresponding index in the gray value sequence. These are the corresponding row / column coordinates; Define loss function : ,in It is the position sequence calculated based on the fitting function during the optimization process. The grayscale value result at the location; Set the optimization objective to Loss convergence; Setting parameters The initial value; Perform nonlinear least squares fitting and iterative optimization. Four parameters, until The optimization stops when the difference from the actual pixel data is minimized, i.e., when the loss value is less than the set threshold or when the set maximum number of iterations is reached; Converged parameters That is, the x-coordinate of the center point of the S-curve; Repeat the above process for all rows / columns to generate a set of sub-pixel point clouds for the four edges. .

4. The real-time correction method based on EPD images according to claim 1, characterized in that, The process of defining a target rectangle based on user-input nominal field parameters, and using a perspective transformation matrix to map the distorted quadrilateral to the target rectangle to complete image correction includes: The vertex coordinates of the target rectangle are determined based on the nominal physical width of the field of view, the nominal physical height of the field of view, and the pixel-to-millimeter coefficient input by the user. Determine the perspective transformation matrix based on the vertex coordinates of the target rectangle and the vertex coordinates of the distorted rectangle. ; Traverse the pixel coordinates in the target image and use a matrix inverse matrix Map the coordinates of the target image back to the source image; A weighted average of the pixels in the source image is generated using a cubic convolution kernel function W(d).

5. The real-time correction method based on EPD images according to claim 4, characterized in that, The perspective transformation matrix for: ； in, For affine components, used to correct image rotation caused by collimator rotation or rotation in the EPD plane, as well as overall scaling caused by source-image distance deviation; The translation component is used to correct the X-axis and Y-axis translation deviations of the field center relative to the imaging center; The perspective component is used to correct the trapezoidal distortion caused by the EPD detector plane not being perpendicular to the beam center axis; These are the normalized parameters.

6. A method for analyzing the correction effect based on EPD images, characterized in that, This includes obtaining a corrected image using the real-time correction method based on EPD images as described in any one of claims 1-5, calculating symmetry based on the corrected image, and calculating confidence based on the difference between the distorted quadrilateral and the corrected image.

7. The method for analyzing the correction effect based on EPD images according to claim 6, characterized in that, The calculation of symmetry based on the corrected image includes: On the corrected image, the horizontal dose distribution curve and the vertical dose distribution curve are extracted along the horizontal and vertical central axes, respectively. Calculate the horizontal and vertical symmetry within the field width n3: Horizontal symmetry is the maximum of the absolute values ​​of the ratio of the difference between symmetrical pixels on the horizontal dose distribution curve to the maximum value on the horizontal dose distribution curve. Vertical symmetry is the maximum of the absolute values ​​of the ratio of the difference between symmetrical pixels on the vertical dose distribution curve to the maximum value on the vertical dose distribution curve. Where n3 is a percentage.

8. The method for analyzing the correction effect based on EPD images according to claim 6, characterized in that, The confidence level calculation based on the difference between the distorted quadrilateral and the corrected image includes: The Y-axis tilt distortion factor is calculated based on the lengths of the top and bottom sides of the distorted quadrilateral. The X-axis tilt distortion factor is calculated based on the left and right lengths of the distorted quadrilateral. Dimensional fidelity is calculated based on the side length of the distorted quadrilateral and the nominal pixel width and nominal pixel height of the corrected image.

9. The method for analyzing the correction effect based on EPD images according to claim 6, characterized in that, It also includes visualizing the symmetry calculation results and confidence calculation results, as well as generating reports.

10. A system for implementing the real-time correction method based on EPD images as described in claim 1, characterized in that, include: The user input module is configured to allow users to input nominal field parameters; The dose image filtering module is configured to automatically filter dose images to obtain cumulative dose images; The field distortion parameter calculation module is configured to locate the four corner points of the distorted field of the electronic field imaging device EPID image based on the cumulative dose image, and obtain the distorted quadrilateral. The perspective transformation correction module is configured to define the target rectangle based on the nominal field parameters input by the user, and use the perspective transformation matrix to map the distorted quadrilateral to the target rectangle, thereby completing the image correction.

11. The system according to claim 10, characterized in that, Also includes: The analysis and confidence calculation module is configured to calculate symmetry based on the corrected image and confidence based on the difference between the distorted quadrilateral and the corrected image; The visualization and reporting module is configured to visualize the symmetry calculation results and confidence calculation results and generate reports.

Images