A method for direct conversion of EPID images to incident fluence for dose reconstruction

A method for direct conversion of EPID images to incident fluence for dose reconstruction 1 2 3 4 first-author 2 3 4 2 3 4 5 1 2 3 4 corresp Corresponding author, liqin.hu@fds.org.cn Corresponding author, liqin.hu@fds.org.cn liqin.hu@fds.org.cn 1 2 3 4 University of Science and Technology of China, Hefei 230027, China University of Science and Technology of China, Hefei 230027, China Key Laboratory of Neutronics and Radiation Safety, Institute of Nuclear Energy Safety Technology, Chinese Academy of Sciences, Hefei 230031, China Key Laboratory of Neutronics and Radiation Safety, Institute of Nuclear Energy Safety Technology, Chinese Academy of Sciences, Hefei 230031, China Anhui Radiotherapy Engineering Technology Research Center, Hefei 230031, China Anhui Radiotherapy Engineering Technology Research Center, Hefei 230031, China Collaborative Innovation Center of Radiation Medicine of Jiangsu Higher Education Institutions, Suzhou 215006, China Collaborative Innovation Center of Radiation Medicine of Jiangsu Higher Education Institutions, Suzhou 215006, China Radiotherapy center, Cancer Hospital, Hefei Institutes of Physical Science, Chinese Academy of Sciences, Hefei 230031, China Radiotherapy center, Cancer Hospital, Hefei Institutes of Physical Science, Chinese Academy of Sciences, Hefei 230031, China A direct incident fluence measurement method based on amorphous silicon electronic portal imaging device (a-Si EPID) has been developed for pretreatment verification of intensity-modulated radiation therapy (IMRT). The EPID-based incident fluence conversion method deconvolves EPID images to the primary response distribution based on measured lateral scatter kernels in the EPID detector using Conjugate Gradient algorithm. The primary response is converted to the incident fluence based on measured fluence conversion matrix which corrects for off-axis position dependence of the a-Si EPID response and the "horn" beam profile caused by flatting filter. To verify feasibility and accuracy of this method, square fields of various sizes and two IMRT plans were delivered. The dose distributions computed based on EPID-derived incident fluence were compared with the measurement data. For all square field sizes except the smallest field (2 cm), the mean dose differences in cross-line dose profiles were within 1% excluding the penumbra region, and gamma passing percentages with a 2%/2 mm criterion were about 99%. For two IMRT plans, the least gamma passing percentage for all eight IMRT fields was 98.14% with 2%/3 mm criteria. It can be concluded that our direct EPID-based incident fluence conversion method is accurate and capable of being applied to pretreatment dose verification in clinical routines. Amorphous-silicon portal imaging devices Dose reconstruction Dose verification I. INTRODUCTION 1 Intensity modulated radiotherapy (IMRT) is a powerful technique to achieve a better dose conformity to the tumor volume and reduce the dose to normal tissues as much as possible. This is only achieved if the planned fluence is delivered accurately at the treatment unit. However, under-dose of the tumor or an over-dose of the organs at risk would be caused by the following factors [ 1 1 ]: (a) the calculated leaf sequence does not result in the fluence pattern used by the treatment plan system for dose calculation; (b) the treatment plan is not correctly transferred to the accelerator; and (c) the treatment equipment is not functioning correctly. Pretreatment verification of the fluence delivered is essential to ensure the accuracy of IMRT treatment. Nevertheless, no detector can be used so far to directly measure the intensity fluence for radiotherapy. Electronic portal imaging device (EPID), an X-ray imaging detector attached to the accelerator, has gone on to replace traditional dosimetry (such as ion chamber, film, etc.) as pretreatment dose verification tool, due to the advantages of fast image acquisition, high resolution, digital format and good dosimetric characteristics [ 2 2 - 5 5 2-5 ]. One way to use the EPID for pretreatment verification is to compare measured image (raw or converted to dose-to-water) with predicted EPID image [ 6 6 - 8 8 6-8 ] or dose-to-water at the level of the EPID [ 9 9 , 10 10 ]. But its limitation is that it cannot directly interpret the discrepancy between the planned doses with the delivered dose inside patient. A more beneficial way to use EPID is to compare 3D dose reconstructed with plan calculated inside patient, provided that the EPID image can be converted to the delivered energy fluence used as input for dose reconstruction. Warkentin et al . [ 11 11 ] deconvoluted the EPID raw image to 2D primary fluence distribution using scatter kernels generated by Monte Carlo model of EPID. Van Elmpt et al . [ 12 12 - 14 14 12-14 ] extracted the energy fluence from a measured dose-to-water with EPID by a back-projection model. Zijtveld et al . [ 15 15 ] used an iterative algorithm to obtain the incident fluence from portal dose image. However, modelling the EPID by Monte Carlo accurately is difficult, as detailed structure of EPID is usually not available. In Refs. [ 12 12 - 15 15 12-15 ], for obtaining the fluence, the EPID images should be converted to portal dose distribution by EPID dose calibration model, a complicate job with easy dosimetry errors. In this paper, a method is proposed to directly convert the EPID image to the incident fluence. It has an effective correction procedure for the field size and off-axis position dependence of the EPID response, based on measured scatter kernels and fluence conversion matrix, which are derived from a series of measurements. EPID grey images are deconvolved to the primary response distributions based on scatter kernels using Conjugate Gradient algorithm. The primary responses are converted to the incident fluence. To validate accuracy of EPID-derived incident fluence, the dose computed based on the measured fluence is compared to measured dose for a large range of conditions. II. MATERIALS AND METHODS 1 A. Electronic portal imaging device (EPID) 2 The EPID is of a PerkinElmer 1640L amporphous silicon type EPID (PerkinElmer Optoelectronics, Wiesbaden, Germany) mounted on a Elekta Axesse linac (Elekta AB, Stockholm, Sweden). On top of the sensitive layer of the EPID, a 3 mm copper slab was mounted to remove low energy photons, which cause over-response of EPID [ 13 13 , 16 16 ]. All images were processed to correct for pixel-to-pixel sensitivity variations and dark current. In addition, image lag and ghosting artifacts affecting amorphous silicon EPIDs [ 4 4 , 17 17 ], were corrected by using a ghosting correction factor G , similar to the procedure described by Nijsten [ 18 18 ]. B. EPID-based incident fluence conversion model 2 1. Determining the scatter kernels 3 The response I ( x , y ) measured at a point ( x , y ) in the EPID consists of two parts: the primary response I P ( x , y ) of radiation reaching the EPID plane directly and the composite scatter response I S ( x , y ) from lateral scatter within the EPID and accelerator head. <math display="block" id="M1"><mrow><mi>I</mi><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mo>=</mo><msub><mi>I</mi><mtext>P</mtext></msub><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mo>+</mo><msub><mi>I</mi><mtext>S</mtext></msub><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo></mrow></math> In our method, the scatter component is modeled as convolution of the primary response with a scatter kernel K ( r ( x-x’ , y-y’ ) , d ): <math display="block" id="M2"><mrow><msub><mi>I</mi><mtext>S</mtext></msub><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mo>=</mo><mstyle displaystyle="true"><mrow><munder><mo>∬</mo><mrow><msup><mi>x</mi><mo>′</mo></msup><mn>,</mn><msup><mi>y</mi><mo>′</mo></msup><mo>∈</mo><mtext>Field</mtext></mrow></munder></mrow></mstyle><msub><mi>I</mi><mtext>P</mtext></msub><mo stretchy="false">(</mo><msup><mi>x</mi><mo>′</mo></msup><mn>,</mn><msup><mi>y</mi><mo>′</mo></msup><mo stretchy="false">)</mo><mo>⋅</mo><mi>K</mi><mo stretchy="false">(</mo><msub><mi>r</mi><mrow><mo stretchy="false">(</mo><mi>x</mi><mo>−</mo><msup><mi>x</mi><mo>′</mo></msup><mn>,</mn><mi>y</mi><mo>−</mo><msup><mi>y</mi><mo>′</mo></msup><mo stretchy="false">)</mo></mrow></msub><mn>,</mn><mi>d</mi><mo stretchy="false">)</mo><mo>⋅</mo><mtext>d</mtext><msup><mi>x</mi><mo>′</mo></msup><mtext>d</mtext><msup><mi>y</mi><mo>′</mo></msup><mn>,</mn></mrow></math> where, I P ( x’ , y’ ) is the primary response at point ( x’ , y’ ), r ( x-x’ , y-y’ ) =[( x-x’ ) 2 +( y-y’ ) 2 ] 1/2 is the distance of point ( x’ , y’ ) to the measured point ( x , y ), K ( r ( x-x’ , y-y’ ) , d ) is related to the off-axis distance d (beam-softening) due to energy dependency of EPID response. In order to derive the scatter kernel, the primary response directly reaching to the EPID is calculated first. Since the scatter component I S ( x , y ) depends on field size, the primary component I P ( x , y ) has to be field size independent. This fact can be used to estimate the primary contribution as follows. The total on-axis response of the EPID I ( x , y ) is experimentally determined as a function of field size f s , by irradiating the EPID with square fields of different sizes. For the limit of zero field size, the scatter component is negligible. Thus, the primary response I P ( x , y ) can be assessed by extrapolating the total response I ( x , y ) to a field size f s of 0 cm×0 cm: <math display="block" id="M3"><mrow><msub><mi>I</mi><mtext>P</mtext></msub><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mo>≈</mo><munder><mstyle mathsize="140%" displaystyle="true"><mrow><mi>lim</mi></mrow></mstyle><mrow><msub><mi>f</mi><mi>s</mi></msub><mo>→</mo><mn>0</mn></mrow></munder><mi>I</mi><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mn>.</mn></mrow></math> From Eqs. (1) and (2), we have: <math display="block" id="M4"><mrow><mtable columnalign="left"><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow><msub><mrow><mrow><mo>(</mo><mrow><mfrac><mrow><mi>I</mi><mo stretchy="false">(</mo><mi>x</mi><mo>=</mo><mn>0,</mn><mi>y</mi><mo>=</mo><mn>0</mn><mo stretchy="false">)</mo><mo>−</mo><msub><mi>I</mi><mtext>P</mtext></msub><mo stretchy="false">(</mo><mi>x</mi><mo>=</mo><mn>0,</mn><mi>y</mi><mo>=</mo><mn>0</mn><mo stretchy="false">)</mo></mrow><mrow><msub><mi>I</mi><mtext>P</mtext></msub><mo stretchy="false">(</mo><mi>x</mi><mo>=</mo><mn>0,</mn><mi>y</mi><mo>=</mo><mn>0</mn><mo stretchy="false">)</mo></mrow></mfrac></mrow><mo>)</mo></mrow></mrow><mrow><mtext>norm</mtext></mrow></msub></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow><mo>=</mo><mstyle displaystyle="true"><mrow><munder><mo>∬</mo><mrow><mi>x</mi><mn>,</mn><mi>y</mi><mo>∈</mo><mtext>Field</mtext></mrow></munder></mrow></mstyle><mfrac><mrow><msub><mi>I</mi><mtext>P</mtext></msub><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo></mrow><mrow><msub><mi>I</mi><mtext>P</mtext></msub><mo stretchy="false">(</mo><mi>x</mi><mo>=</mo><mn>0,</mn><mi>y</mi><mo>=</mo><mn>0</mn><mo stretchy="false">)</mo></mrow></mfrac><mo>⋅</mo><mi>K</mi><mo stretchy="false">(</mo><msub><mi>r</mi><mrow><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo></mrow></msub><mn>,</mn><mi>d</mi><mo stretchy="false">)</mo><mo>⋅</mo><mtext>d</mtext><mi>x</mi><mtext>d</mtext><mi>y</mi></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow><mo>=</mo><mstyle displaystyle="true"><mrow><munder><mo>∬</mo><mrow><mi>x</mi><mn>,</mn><mi>y</mi><mo>∈</mo><mtext>Field</mtext></mrow></munder></mrow></mstyle><mi>f</mi><mo stretchy="false">(</mo><msub><mi>r</mi><mrow><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo></mrow></msub><mo stretchy="false">)</mo><mo>⋅</mo><mi>K</mi><mo stretchy="false">(</mo><msub><mi>r</mi><mrow><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo></mrow></msub><mn>,</mn><mi>d</mi><mo stretchy="false">)</mo><mo>⋅</mo><mtext>d</mtext><mi>x</mi><mtext>d</mtext><mi>y</mi><mn>,</mn></mrow></mtd></mtr></mtable></mrow></math> <math display="block" id="M5"><mrow><mi>K</mi><mo stretchy="false">(</mo><mi>r</mi><mn>,</mn><mi>d</mi><mo stretchy="false">)</mo><mo>=</mo><mrow><mo>{</mo><mrow><mtable columnalign="left"><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow><mstyle displaystyle="true"><munderover><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>0</mn></mrow><mi>N</mi></munderover></mstyle><msub><mi>a</mi><mi>i</mi></msub><mo>⋅</mo><mi>exp</mi><mrow><mo>(</mo><mrow><mo>−</mo><msup><mrow><mrow><mo>(</mo><mrow><mfrac><mi>r</mi><mrow><msub><mi>b</mi><mi>i</mi></msub></mrow></mfrac></mrow><mo>)</mo></mrow></mrow><mn>2</mn></msup></mrow><mo>)</mo></mrow><mtext> </mtext><mtext>for </mtext><mi>r</mi><menclose notation="updiagonalstrike"><mo>=</mo></menclose><mn>0</mn></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow><mn>1</mn><mtext> </mtext><mtext>for </mtext><mi>r</mi><mo>=</mo><mn>0</mn></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd><mtd columnalign="left"><mrow/></mtd></mtr></mtable></mrow></mrow><mn>,</mn></mrow></math> where the point ( x =0, y =0) represent the center of the beam field, f ( r ( x , y ) ) is the off-axis primary response profile, and a i and b i are the kernel parameters. We introduced radial symmetric N Gaussian function to describe the scatter kernel. To accurately derive a i and b i , the scatter contributions for a set of square fields of different sizes are estimated by subtracting the primary response computed using Eq. (3). The scatter data represented as a function of field size are fitted by adjusting the kernel parameters. 2. Primary response extraction 3 Assuming that the incident intensity at the source-axis distance (SAD) of 100 cm is composed of pencil beams of F x × F y finite size, the response for pencil beam ( i , j ) reaching directly to the EPID is the primary response I P ( i , j ). We can derive that the total response value of the point ( x , y ) in the EPID is the sum of the primary response I P ( i , j ) and the scatter component from other pencil beams: <math display="block" id="M6"><mrow><mrow><mo>{</mo><mrow><mtable columnalign="left"><mtr columnalign="left"><mtd columnalign="left"><mrow><msub><mi>I</mi><mi>C</mi></msub><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mo>=</mo><mstyle displaystyle="true"><munderover><mo>∑</mo><mrow><mi>i</mi><mo>=</mo><mn>0</mn></mrow><mrow><msub><mi>F</mi><mi>x</mi></msub></mrow></munderover></mstyle><mstyle displaystyle="true"><munderover><mo>∑</mo><mrow><mi>j</mi><mo>=</mo><mn>0</mn></mrow><mrow><msub><mi>F</mi><mi>y</mi></msub></mrow></munderover></mstyle><msub><mi>I</mi><mi>P</mi></msub><mo stretchy="false">(</mo><mi>i</mi><mn>,</mn><mi>j</mi><mo stretchy="false">)</mo><mo>⋅</mo><mi>K</mi><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo>;</mo><mo>{</mo><mi>i</mi><mn>,</mn><mi>j</mi><mo>}</mo><mo stretchy="false">)</mo></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow><mi>σ</mi><mo>=</mo><msqrt><mrow><mfrac><mn>1</mn><mrow><msub><mi>E</mi><mi>x</mi></msub><mo>⋅</mo><msub><mi>E</mi><mi>y</mi></msub></mrow></mfrac><mo>⋅</mo><mstyle displaystyle="true"><munderover><mo>∑</mo><mrow><mi>x</mi><mo>=</mo><mn>0</mn></mrow><mrow><msub><mi>E</mi><mi>x</mi></msub></mrow></munderover></mstyle><mstyle displaystyle="true"><munderover><mo>∑</mo><mrow><mi>y</mi><mo>=</mo><mn>0</mn></mrow><mrow><msub><mi>E</mi><mi>y</mi></msub></mrow></munderover></mstyle><msup><mrow><mo stretchy="false">(</mo><msub><mi>I</mi><mi>C</mi></msub><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mo>−</mo><msub><mi>I</mi><mi>M</mi></msub><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mo stretchy="false">)</mo></mrow><mn>2</mn></msup></mrow></msqrt></mrow></mtd></mtr><mtr columnalign="left"><mtd columnalign="left"><mrow/></mtd></mtr></mtable></mrow></mrow><mn>,</mn></mrow></math> where I C ( x , y ) and I M( x , y ) are the calculated and measured response values at a point ( x , y ) of the EPID, respectively. E x × E y is the size of EPID image used to reconstruct the primary response I P ( i , j ), and K ( x , y ,{ i , j }) is the scatter contribution for point ( x , y ) in the EPID from pencil beam ( i , j ). Equation (6) is a nonlinear equation with F x × F y variances, so the process of solving it becomes an unconstrained optimization problem. The convergence criterion is based on the difference between the calculated and measured response values on the EPID. An in-house developed Conjugate Gradient algorithm [ 19 19 , 20 20 ] was used to fit all the pencil beam primary response I P ( i , j ) so as to minimize error σ between the calculated and measured values for all the E x × E y pixel points pooled together. This algorithm is efficient for solving large-scale optimization problems and has been applied in ARTS [ 21 21 , 22 22 ] (Accurate/Advanced Radiation Therapy System) developed by the FDS Team, an interdisciplinary team devoting to the research and development of nuclear reactor physics [ 23 23 ], nuclear reactor materials [ 24 24 - 26 26 24-26 ], nuclear reactor engineering [ 27 27 - 29 29 27-29 ], numerical simulation and visualization [ 30 30 ], medical physics and environmental protection, etc. 3. Incident fluence conversion 3 In order to convert the primary response values I P ( i , j ) derived from acquired EPID image to the incident fluence Φ ( x , y ), we have to determine the relationship between the primary response values of the EPID and the number of monitor units (MU), defined as a unit representing the linac output of the accelerator. In our work, the fluence conversion matrix covering the whole detector was established based on the measurement data, which accounted for a general off-axis variation in beam fluence profile and variation in detector response due to changes in off-axis energy spectrum of the beam. We used a CC13 ionization chamber with a Φ 3.76 cm mini-phantom as build-up cap to measure the in-air off-axis radio (OAR) profile at a source-axis distance (SAD) of 100 cm by scanning (IBA Dosimetry GmbH, Schwarzenbruck, Germany) the diagonals of a 25×25 cm 2 reference field. OAR from each half of both diagonal beam profiles was used to calculate an average OAR ( r ( x , y ) ). Under the same field size, EPID images were acquired for the calibrated number of monitor units 100 MU. Then, the primary response value of the EPID for each pencil beam <math id="M7"><mrow><msub><mi>I</mi><mrow><msub><mtext>P</mtext><mrow><mtext>cali</mtext></mrow></msub></mrow></msub><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo></mrow></math> could be reconstructed using Eq. (6). The fluence conversion matrix CU ( x , y ) is determined by dividing the primary response values and the number of monitor units corrected by OAR. <math display="block" id="M8"><mrow><mi>C</mi><mi>U</mi><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mo>=</mo><msub><mi>I</mi><mrow><msub><mtext>P</mtext><mrow><mtext>cali</mtext></mrow></msub></mrow></msub><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mo>/</mo><mo stretchy="false">[</mo><msub><mi>Φ</mi><mrow><mtext>cali</mtext></mrow></msub><mo>⋅</mo><mi>O</mi><mi>A</mi><mi>R</mi><mo stretchy="false">(</mo><msub><mi>r</mi><mrow><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo></mrow></msub><mo stretchy="false">)</mo><mo stretchy="false">]</mo></mrow></math> Finally, the incident fluence Φ ( x , y can be calculated for any field size and shape based on the primary response extracted by Step B and the fluence conversion matrix established by Eq. (7). <math display="block" id="M9"><mrow><mi>Φ</mi><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mo>=</mo><msub><mi>I</mi><mtext>P</mtext></msub><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo><mo>/</mo><mi>C</mi><mi>U</mi><mo stretchy="false">(</mo><mi>x</mi><mn>,</mn><mi>y</mi><mo stretchy="false">)</mo></mrow></math> C. Experimental validation of EPID-derived incident fluence 2 Since direct measurement of incident fluence is difficult, it is impossible to estimate the accuracy of EPID-derived incident fluence by comparing EPID-derived fluence with measured fluence. Instead, comparison of dose calculated using EPID-derived incident fluence to measured dose was performed to verify the EPID-derived incident fluence. An accurate in-house developed dose calculation engine based on Monte Carlo Finite-Size Pencil Beam Model [ 31 31 , 32 32 ] was adopted for dose calculation. 1. Open field verification 3 For verifying the sensitivity and accuracy of the EPID-based fluence conversion procedure for different field sizes and numbers of monitor unit (MU), square field sizes of 2 cm×2 cm, 5 cm×5 cm, 10 cm×10 cm, 15 cm×15 cm and 20 cm×20 cm (100 MU) and monitor units of 2, 5, 10, 50, 100 and 150 MU (10 cm×10 cm) were measured. The absolute doses at the field center and the cross-line dose profiles were measured at a depth of 5 cm using the CC13 chamber in the Blue Phantom water phantom system (IBA Dosimetry GmbH, Schwarzenbruck, Germany) with a source-to-surface distance (SSD) of 100 cm. 2. Clinical cases verification 3 For verifying feasibility of the incident fluence measurement procedure for clinical cases, two step-and-shoot IMRT treatment plans were performed. The IMRT-1 treatment plan was a four-field lung plan with a total of 96 segments. The IMRT-2 treatment plan was a four-field head and neck plan with a total of 80 segments. Two-dimensional (2D) dose distribution was measured at 5 cm depth in a solid water phantom using 2D-ARRAY seven29 dosimetry (PTW-Freiburg, Germany) with an SSD of 95 cm. The 2D dose distribution for each field was evaluated by comparing with the measured dose distribution. III. RESULT AND DISCUSSION 1 A. Scatter kernels 2 Scatter kernels for off-axis positions of 0, 2, 4, 6 and 8 cm were derived from measurements. For the determination of these scatter kernels, square fields of 2, 5, 10, 15 and 20 cm widths were measured with EPID for each off-axis position. The primary responses for different off-axis positions were calculated using Eq. (3). The coefficients a i and b i of the scatter kernel were determined by a fitting procedure using a non-linear least squares method implemented in Matlab (Matlab 7.6.0, The Mathworks, Natick, USA) based on Eqs. (4) and (5). The fitting procedure demonstrated that the best kernel fits for each off-axis position could be obtained if three Gaussian distributions ( N =3) were used to describe the kernels. The scatter kernels for different off-axis position are shown in Fig. 1 Fig. 1 . The fitted parameters of the scatter kernels are given in Table 1 Table 1 . 1 1 Off-axis distance (cm) 1 1 a 1 (cm -2 ) 1 1 a 2 (cm -2 ) 1 1 a 3 (cm -2 ) 1 1 b 1 (cm) 1 1 b 2 (cm) 1 1 b 3 (cm) 1 1 0 1 1 0.02197 1 1 0.003422 1 1 0.4900 1 1 1.601 1 1 4.760 1 1 -0.2509 1 1 2 1 1 0.02184 1 1 0.003370 1 1 0.4886 1 1 1.588 1 1 4.885 1 1 -0.2467 1 1 4 1 1 0.01956 1 1 0.003015 1 1 0.4887 1 1 1.665 1 1 5.272 1 1 -0.2467 1 1 6 1 1 0.02208 1 1 0.003514 1 1 0.4896 1 1 1.591 1 1 4.622 1 1 -0.2503 1 1 8 1 1 0.02352 1 1 0.004250 1 1 0.4899 1 1 1.601 1 1 4.760 1 1 -0.2502 Figure 2 Figure 2 shows that the Gaussian distributions of the scatter kernel for different off-axis positions are close to each other. The small differences are possibly due to the impact of changes in the off-axis energy spectrum and the photons incident direction. The same EPID image was deconvolved to the primary responses using scatter kernels for different off-axis positions in turn. Then, the primary responses were compared and the differences were less than 1%. Thus, the assumption that the field size-dependent scatter kernel is invariant with off-axis position is reasonable, as long as no object is in the beam. B. Experimental validation of EPID-derived incident fluence 2 1. Open field verification 3 The measured and calculated absolute dose profiles in cross-line direction for monitor units of 2, 5, 10, 50, 100 and 150 MU (10 cm×10 cm) and field sizes of 2 cm×2 cm, 5 cm×5 cm, 10 cm×10 cm, 15 cm×15 cm, 20 cm×20 cm (100 MU) are shown in Fig. 2 Fig. 2 . Gamma values are analyzed over the area within 5% isodose line of maximum dose using the criteria of 2%/2 mm. The mean dose differences of 2 and 5 MU are slightly larger than those of monitor units of 10 MU, but they were still within 1.5%, excluding the penumbra region. For the gamma values, 95.5% of 2 MU and 98.7% of 5 MU were less than 1, with the mean gamma values being 0.627 and 0.569, respectively. For the other monitor units of 10 MU, all gamma passing percentages were more than 99% and corresponding mean gamma values were less than 0.4. The larger dose differences for small numbers of monitor units can be possibly explained by the instable output of the linac. Except for 2 cm×2 cm field, the mean dose differences for all field sizes were less than 1% (excluding the penumbra region); gamma passing percentages were about 99% with mean gamma values less than 0.4. The dose error at the center of the field for the 2 cm×2 cm field size was 2.7% and the corresponding gamma passing percentage was 69.2%. The larger dose error for 2 cm×2 cm field size may be caused by a larger uncertainty due to the small field size relative to the ion chamber volume. The mean dose differences become slightly larger if the penumbral regions are included. But the gamma analysis indicates that the calculated doses based on EPID-derived fluence and the measured doses agree quite well within 5% isodose line of maximum dose, including the penumbral regions. This implies that the calculated dose profiles have a steeper dose falloff in the penumbra than the measured dose profiles due to blurring of the measurement by the ionization chamber. 2. Clinical cases verification 3 Two clinical step-and-shoot IMRT plan, consisting of eight fields totally, were measured with EPID and 2D ionization chamber array, respectively. VeriSoft5.1 software (PTW) was used to compare the measured dose distributions with the calculated dose distribution based on EPID-derived fluence. Figure 3 Figure 3 shows the results for one field consisting of 26 segments as an example. The isodose lines ( Fig. 3(a) Fig. 3(a) ) demonstrate very good agreement. Figure 3(b) Figure 3(b) shows the gamma histogram for this field. The gamma distributions with the criteria of 2%/3 cm for all eight IMRT fields were calculated over the area within the 5% isodose line of the maximum dose, which disregard the regions outside the field. The gamma values for all eight IMRT fields are given in Table 2 Table 2 . 1 1 Plan 1 1 Field 1 1 γ mean 1 1 γ max 1 1 P γ < 1 (%) 1 1 IMRT1 1 1 1 1 1 0.385 1 1 1.094 1 1 99.09 1 1 1 1 2 1 1 0.388 1 1 1.112 1 1 99.14 1 1 1 1 3 1 1 0.408 1 1 1.102 1 1 98.78 1 1 1 1 4 1 1 0.396 1 1 0.967 1 1 100.00 1 1 IMRT2 1 1 1 1 1 0.363 1 1 0.952 1 1 100.00 1 1 1 1 2 1 1 0.415 1 1 1.095 1 1 98.48 1 1 1 1 3 1 1 0.396 1 1 0.814 1 1 100.00 1 1 1 1 4 1 1 0.425 1 1 1.132 1 1 98.14 All gamma distributions had a mean value below 0.43. The least gamma passing percentage ( P <1 ) for all IMRT fields was 98.14%. The gamma passing percentages for three fields reached 100%, implying that all the measured points satisfied the chosen gamma criteria. Excellent agreement between computed dose based on EPID-derived incident fluence and measured dose indicates that the incident fluence obtained by our method is accurate. The reconstructed incident fluence can be used as input to TPS, so as to generate a 3D dose distribution using the patient’s CT data. It can provide the cumulative errors in a 3D dose distribution in the patient from all beams in the IMRT plan to help the clinician to assess the potential clinical significance of dosimetry uncertainties. Thus, if the method can be streamlined sufficiently, 3D dose verification may be a clinically useful quality assurance tool that provides information not easily accessible using more conventional 2D verification method. IV. CONCLUSION 1 In this article, we have developed an accurate EPID-based incident fluence conversion procedure to directly convert EPID images to the incident fluence. For this purpose, a conversion model has been investigated that corrects for the field size and the off-axis position dependence of the EPID response based on measured scatter kernels and fluence conversion matrix, without the need of input of beam characteristics and detailed analysis of EPID construction. EPID-derived incident fluences were evaluated for different fields such as open fields and IMRT fields to determine the accuracy and applicability of the EPID-based incident fluence conversion model. Gamma evaluation criteria were satisfied for all measurements. Using our EPID-based incident fluence conversion method, highly accurate incident fluence can be obtained which are promising for pretreatment 3D dose verification. References M Van Zijtveld , M L P Dirkx , H C J De Boer , et al . Dosimetric pre-treatment verification of IMRT using an EPID; clinical experience . Radiother Oncol , 2006 , 81 : 168 - 175 . DOI: 10.1016/j.radonc.2006.09.008 http://doi.org/10.1016/j.radonc.2006.09.008 Dosimetric pre-treatment verification of IMRT using an EPID; clinical experience P B Greer . Off-axis dose response characteristics of an amorphous silicon electronic portal imaging device . Med Phys , 2007 , 34 : 3815 - 3824 . DOI: 10.1118/1.2779944 http://doi.org/10.1118/1.2779944 Off-axis dose response characteristics of an amorphous silicon electronic portal imaging device R J Louwe , L N McDermott , J J Sonke , et al . The long-term stability of amorphous silicon flat panel imaging devices for dosimetry purposes . Med Phys , 2004 , 31 : 2989 - 2995 . DOI: 10.1118/1.1803751 http://doi.org/10.1118/1.1803751 The long-term stability of amorphous silicon flat panel imaging devices for dosimetry purposes L N McDermott , R J W Louwe , J J Sonke , et al . Dose-response and ghosting effects of an amorphous silicon electronic portal imaging device . Med Phys , 2004 , 31 : 285 - 295 . DOI: 10.1118/1.1637969 http://doi.org/10.1118/1.1637969 Dose-response and ghosting effects of an amorphous silicon electronic portal imaging device P Winkler , A Hefner and D Georg . Dose-response characteristics of an amorphous silicon EPID . Med Phys , 2005 , 32 : 3095 - 3105 . DOI: 10.1118/1.2040711 http://doi.org/10.1118/1.2040711 Dose-response characteristics of an amorphous silicon EPID L Parent , A L Fielding , D R Dance , et al . Amorphous silicon EPID calibration for dosimetric applications: comparison of a method based on Monte Carlo prediction of response with existing techniques . Phys Med Biol , 2007 , 52 : 3351 - 3368 . DOI: 10.1088/0031-9155/52/12/003 http://doi.org/10.1088/0031-9155/52/12/003 Amorphous silicon EPID calibration for dosimetric applications: comparison of a method based on Monte Carlo prediction of response with existing techniques L Parent , J Seco , P M Evans , et al . Monte Carlo modelling of a-Si EPID response: The effect of spectral variations with field size and position . Med Phys , 2006 , 33 : 4527 - 4540 . DOI: 10.1118/1.2369465 http://doi.org/10.1118/1.2369465 Monte Carlo modelling of a-Si EPID response: The effect of spectral variations with field size and position J V Siebers , J O Kim , L Ko , et al . Monte Carlo computation of dosimetric amorphous silicon electronic portal images . Med Phys , 2004 , 31 : 2135 - 2146 . DOI: 10.1118/1.1764392 http://doi.org/10.1118/1.1764392 Monte Carlo computation of dosimetric amorphous silicon electronic portal images J Chang , G S Mageras , C S Chui , et al . Relative profile and dose verification of intensity-modulated radiation therapy . Int J Radiat Oncol , 2000 , 47 : 231 - 240 . DOI: 10.1016/S0360-3016(99)00555-6 http://doi.org/10.1016/S0360-3016(99)00555-6 Relative profile and dose verification of intensity-modulated radiation therapy A Van Esch , T Depuydt and D P Huyskens . The use of an a Si-based EPID for routine absolute dosimetric pre-treatment verification of dynamic IMRT fields . Radiother Oncol , 2004 , 71 : 223 - 234 . DOI: 10.1016/j.radonc.2004.02.018 http://doi.org/10.1016/j.radonc.2004.02.018 The use of an a Si-based EPID for routine absolute dosimetric pre-treatment verification of dynamic IMRT fields B Warkentin , S Steciw , S Rathee , et al . Dosimetric IMRT verification with a flat-panel EPID . Med Phys , 2003 , 30 : 3143 - 3155 . DOI: 10.1118/1.1625440 http://doi.org/10.1118/1.1625440 Dosimetric IMRT verification with a flat-panel EPID W J C Van Elmpt , S M J J G Nijsten , R F H Schiffeleers , et al . A Monte Carlo based three-dimensional dose reconstruction method derived from portal dose images . Med Phys , 2006 , 33 : 2426 - 2434 . DOI: 10.1118/1.2207315 http://doi.org/10.1118/1.2207315 A Monte Carlo based three-dimensional dose reconstruction method derived from portal dose images W Van Elmpt , S Petit , D De Ruysscher , et al . 3D dose delivery verification using repeated cone-beam imaging and EPID dosimetry for stereotactic body radiotherapy of non-small cell lung cancer . Radiother Oncol , 2010 , 94 : 188 - 194 . DOI: 10.1016/j.radonc.2009.12.024 http://doi.org/10.1016/j.radonc.2009.12.024 3D dose delivery verification using repeated cone-beam imaging and EPID dosimetry for stereotactic body radiotherapy of non-small cell lung cancer W Van Etmpt , S Nijsten , B Mijnheer , et al . The next step in patient-specific QA: 3D dose verification of conformal and intensity-modulated RT based on EPID dosimetry and Monte Carlo dose calculations . Radiother Oncol , 2008 , 86 : 86 - 92 . DOI: 10.1016/j.radonc.2007.11.007 http://doi.org/10.1016/j.radonc.2007.11.007 The next step in patient-specific QA: 3D dose verification of conformal and intensity-modulated RT based on EPID dosimetry and Monte Carlo dose calculations M Van Zijtveld , M L P Dirkx , H C J De Boer , et al . 3D dose reconstruction for clinical evaluation of IMRT pretreatment verification with an EPID . Radiother Oncol , 2007 , 82 : 201 - 207 . DOI: 10.1016/j.radonc.2006.12.010 http://doi.org/10.1016/j.radonc.2006.12.010 3D dose reconstruction for clinical evaluation of IMRT pretreatment verification with an EPID F Cremers , T Frenzel , C Kausch , et al . Performance of electronic portal imaging devices (EPIDs) used in radiotherapy: Image quality and dose measurements . Med Phys , 2004 , 31 : 985 - 996 . DOI: 10.1118/1.1688212 http://doi.org/10.1118/1.1688212 Performance of electronic portal imaging devices (EPIDs) used in radiotherapy: Image quality and dose measurements L N McDermott , S M J J G Nijsten , J J Sonke , et al . Comparison of ghosting effects for three commercial a-Si EPIDs . Med Phys , 2006 , 33 : 2448 - 2451 . DOI: 10.1118/1.2207318 http://doi.org/10.1118/1.2207318 Comparison of ghosting effects for three commercial a-Si EPIDs S M J J G Nijsten , W J C van Elmpt , M Jacobs , et al . A global calibration model for a Si EPIDs used for transit dosimetry . Med Phys , 2007 , 34 : 3872 - 3884 . DOI: 10.1118/1.2776244 http://doi.org/10.1118/1.2776244 A global calibration model for a Si EPIDs used for transit dosimetry R F Cao , Y C Wu , X Pei , et al . Multi-objective optimization of inverse planning for accurate radiotherapy . Chinese Phys C , 2011 , 35 : 313 - 317 . DOI: 10.1088/1674-1137/35/3/019 http://doi.org/10.1088/1674-1137/35/3/019 Multi-objective optimization of inverse planning for accurate radiotherapy G Li , G Song , Y C Wu , et al . A multi-objective hybrid genetic based optimization for external beam radiation . Plasma Sci Technol , 2006 , 8 : 234 - 236 . A multi-objective hybrid genetic based optimization for external beam radiation S X Tao , Y C Wu , Y X Chen , et al . Repeated positioning in accurate radiotherapy based on virtual net technique and contrary reconstruction scheme . Comput Med Imag Grap , 2006 , 30 : 273 - 278 . DOI: 10.1016/j.compmedimag.2006.05.006 http://doi.org/10.1016/j.compmedimag.2006.05.006 Repeated positioning in accurate radiotherapy based on virtual net technique and contrary reconstruction scheme Y Wu , G Song , C Ruifen , et al . Development of accurate/advanced radiotherapy treatment planning and quality assurance system (ARTS) . Chinese Phys C , 2008 , 32 : 177 - 182 . Development of accurate/advanced radiotherapy treatment planning and quality assurance system (ARTS) Y Wu , Z Xie and U Fischer . A discrete ordinates nodal method for one-dimensional neutron transport calculation in curvilinear geometries . Nucl Sci Eng , 1999 , 133 : 350 - 357 . A discrete ordinates nodal method for one-dimensional neutron transport calculation in curvilinear geometries Q Y Huang , S Gao , Z Q Zhu , et al . Progress in compatibility experiments on lithium-lead with candidate structural materials for fusion in China . Fusion Eng Des , 2009 , 84 : 242 - 246 . DOI: 10.1016/j.fusengdes.2008.12.038 http://doi.org/10.1016/j.fusengdes.2008.12.038 Progress in compatibility experiments on lithium-lead with candidate structural materials for fusion in China Q Huang , C Li , Y Li , et al . Progress in development of China low activation martensitic steel for fusion application . J Nucl Mater , 2007 , 367 : 142 - 146 . DOI: 10.1016/j.jnucmat.2007.03.153 http://doi.org/10.1016/j.jnucmat.2007.03.153 Progress in development of China low activation martensitic steel for fusion application Q Wu , Z Huang , Z Zhu , et al . R&D of dragon series lithium lead loops for material and blanket technology testing . Fusion Sci Technol , 2012 , 62 : 272 - 275 . R&D of dragon series lithium lead loops for material and blanket technology testing Y C Wu . Conceptual design and testing strategy of a dual functional lithium-lead test blanket module in ITER and EAST . Nucl Fusion , 2007 , 47 : 1533 - 1539 . DOI: 10.1088/0029-5515/47/11/015 http://doi.org/10.1088/0029-5515/47/11/015 Conceptual design and testing strategy of a dual functional lithium-lead test blanket module in ITER and EAST Y C Wu , J Q Jiang , M H Wang , et al . A fusion-driven subcritical system concept based on viable technologies . Nucl Fusion , 2011 , 51 : 103036 . DOI: 10.1088/0029-5515/51/10/103036 http://doi.org/10.1088/0029-5515/51/10/103036 A fusion-driven subcritical system concept based on viable technologies Y C Wu , J P Qian and J N Yu . The fusion-driven hybrid system and its material selection . J Nucl Mater , 2002 , 307 : 1629 - 1636 . DOI: 10.1016/s0022-3115(02)01272-2 http://doi.org/10.1016/s0022-3115(02)01272-2 The fusion-driven hybrid system and its material selection Y C Wu . CAD-based interface programs for fusion neutron transport simulation . Fusion Eng Des , 2009 , 84 : 1987 - 1992 . DOI: 10.1016/j.fusengdes.2008.12.041 http://doi.org/10.1016/j.fusengdes.2008.12.041 CAD-based interface programs for fusion neutron transport simulation G Li , H Lin , A D Wu , et al . Realization and comparison of several regression algorithms for electron energy spectrum reconstruction . Chin Phys Lett , 2008 , 25 : 2710 - 2713 . Realization and comparison of several regression algorithms for electron energy spectrum reconstruction H Q Zheng , G Y Sun , G Li , et al . Photon dose calculation method based on Monte Carlo finite-size pencil beam model in accurate radiotherapy . Commun Comput Phys , 2013 , 14 : 1415 - 1422 . DOI: 10.4208/cicp.221212.100413a http://doi.org/10.4208/cicp.221212.100413a Photon dose calculation method based on Monte Carlo finite-size pencil beam model in accurate radiotherapy The authors would like to thank the great help from FDS Team. 10.13538/j.1001-8042/nst.26.050201 Supported by National Natural Science Foundation of China (No. 11305203) and Natural Science Foundation of Anhui Province (No. 1508085QH180) 2014-10-24 2015-03-13 2015-10-20 2015-10-20 2021-04-11 2015 Nuclear Science and Techniques 5 26