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

LOW ENERGY ACCELERATOR, RAY AND APPLICATIONS

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

REN Qiang，

CAO Rui-Fen，

PEI Xi，

LI Bing-Bing，

HU Li-Qin，

WU Yi-Can

Nuclear Science and Techniques

Vol.26, No.5

Article number 050201

Published in print 20 Oct 2015

Available online 20 Oct 2015

DOI：10.13538/j.1001-8042/nst.26.050201

38700

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 devicesDose reconstructionDose verification

I. INTRODUCTION

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]: (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-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-8] or dose-to-water at the level of the EPID [9, 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] 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-14] extracted the energy fluence from a measured dose-to-water with EPID by a back-projection model. Zijtveld et al. [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-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

A. Electronic portal imaging device (EPID)

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, 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, 17], were corrected by using a ghosting correction factor G, similar to the procedure described by Nijsten [18].

B. EPID-based incident fluence conversion model

1. Determining the scatter kernels

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.

I (x, y) = I_{P} (x, y) + I_{S} (x, y)

(1)

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):

I_{S} (x, y) = \iint_{x^{'}, y^{'} \in Field} I_{P} (x^{'}, y^{'}) \cdot K (r_{(x - x^{'}, y - y^{'})}, d) \cdot d x^{'} d y^{'},

(2)

where, I_P(x’,y’) is the primary response at point (x’,y’), r_{(x-x’,y-y’)}=[(x-x’)²+(y-y’)²]^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 fs, 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 fs of 0 cm×0 cm:

I_{P} (x, y) \approx \lim_{f_{s} \to 0} I (x, y) .

(3)

From Eqs. (1) and (2), we have:

\begin{array}{l} {(\frac{I (x = 0, y = 0) - I_{P} (x = 0, y = 0)}{I_{P} (x = 0, y = 0)})}_{norm} \\ = \iint_{x, y \in Field} \frac{I_{P} (x, y)}{I_{P} (x = 0, y = 0)} \cdot K (r_{(x, y)}, d) \cdot d x d y \\ = \iint_{x, y \in Field} f (r_{(x, y)}) \cdot K (r_{(x, y)}, d) \cdot d x d y, \end{array}

(4)

K (r, d) = {\begin{array}{l} \sum_{i = 0}^{N} a_{i} \cdot \exp (- {(\frac{r}{b_{i}})}^{2}) for r = 0 \\ 1 for r = 0 \end{array},

(5)

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 ai and bi are the kernel parameters. We introduced radial symmetric N Gaussian function to describe the scatter kernel. To accurately derive ai and bi, 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

Assuming that the incident intensity at the source-axis distance (SAD) of 100 cm is composed of pencil beams of Fx×Fy 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:

{\begin{array}{l} I_{C} (x, y) = \sum_{i = 0}^{F_{x}} \sum_{j = 0}^{F_{y}} I_{P} (i, j) \cdot K (x, y; {i, j}) \\ σ = \sqrt{\frac{1}{E_{x} \cdot E_{y}} \cdot \sum_{x = 0}^{E_{x}} \sum_{y = 0}^{E_{y}} {(I_{C} (x, y) - I_{M} (x, y))}^{2}} \end{array},

(6)

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. Ex× Ey 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 Fx×Fy 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, 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 Ex×Ey pixel points pooled together. This algorithm is efficient for solving large-scale optimization problems and has been applied in ARTS [21, 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], nuclear reactor materials [24-26], nuclear reactor engineering [27-29], numerical simulation and visualization [30], medical physics and environmental protection, etc.

3. Incident fluence conversion

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² 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 $I_{P_{cali}} (x, y)$ 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.

C U (x, y) = I_{P_{cali}} (x, y) / [Φ_{cali} \cdot O A R (r_{(x, y)})]

(7)

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).

Φ (x, y) = I_{P} (x, y) / C U (x, y)

(8)

C. Experimental validation of EPID-derived incident fluence

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, 32] was adopted for dose calculation.

1. Open field verification

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

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

A. Scatter kernels

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 ai and bi 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. The fitted parameters of the scatter kernels are given in Table 1.

Fig. 1.

Scatter kernels for off-axis position of 0–8 cm.

Fitted parameters of scatter kernels for different off-axis positions

Off-axis distance (cm)	a₁ (cm^-2)	a₂ (cm^-2)	a₃ (cm^-2)	b₁ (cm)	b₂ (cm)	b₃ (cm)
0	0.02197	0.003422	0.4900	1.601	4.760	-0.2509
2	0.02184	0.003370	0.4886	1.588	4.885	-0.2467
4	0.01956	0.003015	0.4887	1.665	5.272	-0.2467
6	0.02208	0.003514	0.4896	1.591	4.622	-0.2503
8	0.02352	0.004250	0.4899	1.601	4.760	-0.2502

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.

Fig. 2.

The absolute dose profiles for different number of monitor units in 10 cm×10 cm (a) and different field sizes in 100 MU (b). Solid lines, measured dose; dotted lines, computed dose based on the EPID-derived fluence.

B. Experimental validation of EPID-derived incident fluence

1. Open field verification

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. 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

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 shows the results for one field consisting of 26 segments as an example. The isodose lines (Fig. 3(a)) demonstrate very good agreement. 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.

Fig. 3.

(Color online) Isodose lines of the measured dose (solid lines) and computed dose (dashed lines) based on the EPID-derived fluence for pretreatment verification of an IMRT field (a) and histograms of the gamma distribution (b).

Analysis of gamma-values for all IMRT fields using VeriSoft software (PTW)

Plan	Field	γ_mean	γ_max	P_γ<1(%)
IMRT1	1	0.385	1.094	99.09
	2	0.388	1.112	99.14
	3	0.408	1.102	98.78
	4	0.396	0.967	100.00
IMRT2	1	0.363	0.952	100.00
	2	0.415	1.095	98.48
	3	0.396	0.814	100.00
	4	0.425	1.132	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

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

[1]

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