GB2360684A - Analysis of portal image for radiotherapy treatment - Google Patents

Abstract

The analysis involves a check on the apparatus setup and a check on the positioning of the patient. To check the apparatus setup a reference radiation field image is compared with a sample radiation field image to determine a difference image. At least the sample image is filtered to remove low frequency variations. A region of interest in the reference and sample images is selected and correlated with each other to determine relative displacements of the reference and sample image. The patient positioning is checked by determining relative displacements between reference image and sample image. The invention can be used for positioning of patient for radiation treatment of cancer to determine if radiation beam is correctly set up and if patient is positioned correctly, analysis of portal image is derived during first treatment session or using digitally reconstructed radiograph.

Description

2360684 1 ANALYSIS OF RADIOGRAPHIC IMAGES The present invention relates to

the analysis of radiographic images.

Linear accelerators have been in use for a significant time for the treatment of severe medical conditions such as cancer. The principle of operation of these devices is essentially that electromagnetic radiation of sufficiently high energy (greater than 1MeV) causes the death andlor destruction of living cells in its path. If therefore a beam of high energy photons or electrons is directed to a cancerous area, the cancer will be destroyed.

It is clear that this treatment has the potential to cause harm to the patient if the radiation is incident on healthy tissue. It is therefore common practice to collimate the beam in order to limit its spatial extent and therefore limit the amount of tissue, healthy or otherwise, upon which it is incident. Collimators to this effect are known, and include collimators of variable size and shape which are capable of selectively blocking the radiation beam at its edges thereby to produce a beam whose shape corresponds to that of the tumour.

In this situation, it is clearly necessary to ensure that the positioning of the patient is accurate. This is generally done by shining one or more collimated light sources on the patient, such as low energy lasers. When the patient is first treated, ink marks are placed on the patient's skin and these 2 are subsequently used to align the patient relative to the lasers. This is able to provide a positioning accuracy of approximately 2-5mm.

As a result, the shaped radiation beam needs to be formed to the approximate outline shape of the tumour plus a margin in each direction to allow for positioning error. This results in a significant volume of tissue around the tumour which is needlessly irradiated. If the patient could be positioned more accurately beneath the beam, then this margin around the tumour could be reduced thereby limiting the damage caused to the patient's healthy tissue.

It is known to prepare a "Portal image" derived from the radiation transmitted through the patient and incident on a photographic plate or other detector. However, the absorption co-efficients of bone, tissue etc. at the energies normally used for radiotherapy are very similar and therefore a Portal image shows very little detail.

The present invention provides a method of electronically analysing the Portal image which is able to distil from the initially indistinct image information which is relevant to the treatment and to the positioning of the patient. In a first general stage, the method determines whether or not the radiation beam is correctly set up. This first stage is described and claimed in GB 2331440, from which this application is divided. In the second general stage, according to the present invention, the method deterffiines whether or not the patient is positioned correctly in the radiation beam.

Thus, the present invention proposes a method of distilling from the Portal images sufficient information regarding the patient's internal structure to allow determination of the patient positioning within the radiation field.

It therefore provides a method of analysing a radiographic image comprising the following steps:

(a) preparing a reference image; (b) preparing a sample image; (c) filtering at least the sample image to remove low frequency variations; (d) selecting a region of interest (R01) within the reference and sample images; 3 (e) correlating the ROI of the reference image with the ROI of the sample image; (f) determining the relative displacements of the reference image and the sample image.

The region of interest (R01) may be a selected field of the image.

Alternatively, it can be the entire image, although it is particularly preferred if the edge of the image is avoided. This is because the edge generally includes details of the radiation field edge, which tends to be either very bright or very dark. Thus, such areas of the image tend to predominate in the correlation step. It is therefore preferred if a border is eliminated by the steps of shrinking the radiation field image generated in the first gen ' eral stage. This shrunken image can, optionally, be thresholded and multiplied with the sample andlor reference image of the second general stage.

Shrinking a bit map image can easily be carried out by setting the pixel to a dark state unless all four neighbours are bright. To shrink an image by more than one pixel, the process can either be repeated, or can be carried out by examining pixels a set distance away in each direction. In the present invention, it is preferred if that distance is between 15 and 30 pixels.

It is naturally preferred if the invention is carried out in combination with the first general stage, and the output result is the calculated difference between the radiation field displacement and the image displacement. This will reveal the displacement of the patient relative to the radiation field, which is of course the value of interest.

According to the present invention, the reference image can be either a Portal image from the first treatment, or a pre-treatment test image, or a digitally reconstructed radiograph.

The sample image can be a Portal image taken during a radiotherapy treatment.

It is particularly preferred if the sample image is captured electronically and analysed during the therapeutic radiotherapy treatment. Such treatments usually last one minute or two minutes, and sufficient data is usually available to produce a sample image within up to 15 seconds. It is 4 preferred if data is captured during the first 10 seconds or less, particularly preferred if it is captured during the first 5 seconds or less.

It would be advisable to end the radiotherapy session if inspection according to either stage reveals a significant error. This will allow the patient to be repositioned more accurately, preventing unnecessary irradiation of healthy tissue.

The present invention also relates to a radiotherapy treatment apparatus comprising means for digitally capturing a Portal image and means for analysis of that image, adapted to carry out the steps of any image analysis method set out above.

An embodiment of the present invention will now be described by way of example, with reference to the accompanying Figures; in which Figure 1 shows a histogram of a typical Portal image; Figure 2 shows an image after thresholding; Figure 3 shows a difference image; Figure 4 shows a sample image after suppression of high frequency variations (noise); Figure 5 shows the result of eliminating all but slowly varying greylevels from Figure 4; Figure 6 illustrates the Bandpass filter of the present invention; Figure 7 shows Figure 5 after removal of the border; Figure 8 shows the image of Figure 7 after inversion; Figure 9 shows the correlation co-efficient of two grey-level images prior to analysis according to the present invention; Figure 10 shows the correlation co-efficient of two images after analysis according to the present invention; Figure 11 shows a digitally reconstructed radiograph; Figure 12 shows a difference image of two Portal images prior to anatomical registration; Figure 13 shows a difference image of two Portal images after anatomical registration; Figure 14 shows a difference image of a DIRR and a Portal image p'riQr to anatomical registration; and Figure 15 shows a difference image of a DRIR and a Portal image after anatomical registration.

To assist in a complete understanding of the present invention, the following description contains details of the complete image analysis method, including both the first general stage (which is the subject of the application referred to. above from which this application is divided) and also the second general stage which is the subject of the present invention.

An analysis method according to the present invention will noVIV be described. The method is divided into the following steps:

1. Elimination of detector-related point errors in the Portal image. 2. Segmentation of the radiation field and calculation of the relative translation, rotation and scaling of the radiation fields including verification of geometrical similarity.

3. Automatic Contrast Enhancement 4. Noise suppression by a Gaussian filter. 5. Elimination of the slowly varying grey-levels.

6 6. Elimination of the border of the radiation field. 7. Determination of the translation and rotation of the anatomy by calculating the correlation between the two pre-processed images with or without an ROL In this discussion, the reference image is a Portal image derived during the first treatment session. It is however possible according to the invention to replace this image with a digitally reconstructed radiograph as illustrated in Figure 8. Such radiographs are known per se, together with method support their construction. The principle of the invention is the same regardless of which type of images used to create the reference image, and therefore no further discussion will be given. Other reference images can be used, such as a diagnostic X-ray image from a treatment simulation, or a hard drawing of bone structures of the patient. A DIRR or X-ray image has the advantage of a higher contrast as compared to a Portal image.

1. Elimination of Point Errors Portal images contain both noise and point errors, ie very dark or very bright pixels with grey-levels which are markedly different compared with their neighbours. These arise when the camera receives a direct hit from a scattered MeV photon. They are not therefore an informative part of the image. Such errors should ideally be eliminated, although this step may not be necessary depending on the quality of the image capture apparatus. A Gaussian filter will not eliminate such an error as it serves only to reduce the peak height of a point error and smear the error out amongst neighbouring pixels. Other linear filters generally have a similar effect. Elimination of point errors can in general only be achieved with a non- linear filter such as a median filter. A 3 x 3 median filter will in the present case achieve the desired result. In a 3 x 3 median filter, the value of the centre pixel is replaced with the median value of its eight neighbours. Thus, individual point errors are eliminated completely with little effect on other types of grey-level variation.

7 2. Segmentation and Matching of the Radiation Field

The histograms of the two median filtered grey-level images from Section 1 are calculated. In the case of Portal images, the histogram will generally be bimodal, one maximum belonging to the radiation field and the other belonging to the part of the image which lies in the shadow of the collimator. A typical histogram is shown at Figure 1. The mean position between the two maxima is used as a threshold for segmentation of the radiation field. Thus, each pixel is essentially thresholded according to whether its individual grey-level is above or below the threshold. This generates a binary image of the radiation field in which grey-levels are either 1 or 0. Thus, the binary image bfij) of the radiation field is calculated by setting pixels with a grey-level ftij) larger than the threshold to 1 and pixels with a grey-level f(ij) smaller than the threshold to 0; b (i, j) 1, if f (1,j) >- Threshold { 0, if f (1, j) < Threshold Such a binary image is illustrated at Figure 2, and essentially shows the collimator shape.

After this has been done to both the sample and the reference image, the relative translation. rotation and scaling of two radiation fields is determined straightforwardly by calculating the moments of their binary images. The first order moments of the binary image bMJ) deliver the centre of gravity (xy) of the radiation field:

1.b (1, j) j-b (1, j) X - andy- The orientation and size of the radiation field can be inferred from the second order moment matrix Mk, which is defined as follows:

(I- X)k. (j -Y) b (1, j) Mki E b (i, j) with k,l=0,1,2 and k+l=2 1,j 8 The matrix M02 M11 M11 M20 is diagonalized and its eigenvalues and corresponding eigenvectors are determined. One gets for the larger eigenvalue 2 -M21 +4 1 M20+M02+ 2M02 2 The angle 0 of its eigenvector with the x-axis is:

M -M + _M20) 2 +4 - M2 e=arctan 02 20 (M02 11 2 M,, The relative translation (Ax,Ay), relative rotation AE) and scaling a of the radiation fields are then:

AX=x portal - X referencer AY= Yportal-yreference AE)=8 portal _e reference X Portal Xreference Usually, large values for the relative translation, rotation and scaling indicated that the two radiation fields are of different shape whereas small values indicate that the radiation field has been chosen correctly. Nevertheless, a large value, e.g. for the translation does not necessarily mean that the radiation field is wrong.

9 To detect errors in the radiation field, one of the binary images is back transformed (inverse translation, inverse rotation and inverse scaling) and then subtracted from the other in order to verify the geometrical similarity of the two images. If both radiation fields are almost identical as applied to the patient then nearly all pixels of the difference image thus generated will be zero, as shown in Figure 3. This last step is valuable because a simple translation of the radiation field will be meaningless, for example showing a translation of the camera between capturing the two images. This could, for example, be caused by removal and refitting of the camera between the two treatments. In that case a translation of the radiation field will be apparent despite the patient being in the correct position relative to the beam. However, demonstrating the geometric similarity of the collimated sha'pqs demonstrates that the apparatus has been configured correctly and is applying the correct radiation field to the patient. Errors at this level are easily made but very significant in terms of patient outcome.

3. Automatic Contrast Enhancement As mentioned above, the histogram of the median-filtered grey-level image has two peaks. The right peak (as shown in Figure 1) contains the anatomical information. To enhance the anatomical structures of interest, an automatic contrast enhancement is performed: all grey-levels on the left side of that peak are set to zero, all grey-levels on the right side are set to 255 and in between a linear grey-level transformation is applied. The ripsult is an image where bony structures are much more visible to a human observer.

4. Noise Suppression A Gaussian filter is applied for noise suppression of the grey-level images demonstrated at Section 1. This leads to an improvement of the signal to noise ratio through the suppression of high frequencies, which usually constitute noise. Gaussian filters of the size between W and 15x '15 pixels generally give good results.

5. Elimination of the Slowly Varying Grey-levels This step is essentially a calculation of the rate of change of pixels, thereby eliminating slow variations. This step leads from a grey-level based correlation method to a structure based correlation method which gives substantially more robust results. If slowly varying grey-levels are allowed to remain in the image, these have a significant influence on the correlation. However, after processing the image from Section 4 above (Figure 4) with a 1 5xl 5 block filter, an image is generated which contains only information regarding slowly varying grey-levels. If this image is subtracted from the image produced at Section 4 (Figure 4), an image is produced where edges are relatively enhanced (Figure 5) and slowly varying grey-levels are eliminated. This is in essence a high pass filter. Thus, in combination with Section 33, a high pass and a low pass filter are applied to the image. This order could clearly be reversed with little detrimental effect on the image.

It is of course the positions of the edges embedded in the image which is necessary for image registration; air bubbles are very distinct in nonstructure based images and tend to influence correlation; by the above process they are reduced from an area-like structure to simply bare border lines which are then almost negligible in correlation. Thus, a matching is possible despite the presence of inevitable air bubbles.

This can be achieved by convoluting the contrast enhanced im"iges with a Bandpass filter to extract bony structures as much as possible. Noise reduction is done by the suppression of high frequencies. The suppression of low frequencies removes the slowly varying background intensities and reduces air bubbles from an area-like structure to their borderlines which are almost negligible in the correlation. The Bandpass filter is a combination of a blockfilter b(xy) filtering with a square of constant pixel values) and a Gaussian f ilter g(x, y). The f iltered image p(x, y) is given by p (x, y) = gfx, y) (1(x, y)-b(x, y)).f(x, y) where the symbol means convolution and Uxy) is the identity operation. After testing with hundreds of image pairs from different body regions, it was found that 1 5xl 5 convolution kernels for b(xy) and g(xy) are good choices for the Bandpass filter (since the pixel size of the Portal images we used is 0.5mm, this results in a mean-frequency of roughly 0.3mm-'. Because of the relatively large kernel sizes the convolution was carried out in Fourier domain by a Fast Fourier Transformation (FFT). Figure 6 illustrates the Bandpass filter.

6. Elimination of the Border The elimination of the border from the image is easiest if the radiation field image (Section 2) is employed. This radiation field image can be shrunk by examining each point in the radiation field to determine whether any one neighbouring point with the distance a (ie (x + a,y), (x-a,y), (x,y + a), (x,y-a)) equals zero. If this is the case, the first point is set to zero in the ry ia,-k image. If all four neighbours have the value 1, the point retains its original value. The result of this is a shrunken radiation field image, which after multiplication leads to the image in Figure 5. This is essentially the same but without the border. This border, consisting of a prominent black and/or white area would otherwise contribute significantly to the correlation.

The distance A can be approximately 30 pixels when a 15x '15 block filter is used, ie about twice the size of the block filter. Suitably, the ratio is between 0.5 and 3 times the size of the block filter.

7. Determination of the Translation and Rotation of the Anatomy In a reference image from Section 6 (Figure 8), one chooses a region of interest (R01) which contains the anatomical structu res of interest and does not contain artifacts. It is often possible to use the entire reference image for determination, but this is not necessary and in particular cases the skill of the radiographer can be deploy in selecting an ROI which gives the best correlation. A major advantage of the present invention is that reliable results are possible based on the entire image; this permits automatic of the comparison and to off-line matching. The similarity of the ROI with a shifted section of the second Portal image from Section 5 is determined by calculating the normalised cross correlation (NCC) as a suitable similarity measure. The maximum of the NCC thus determines the translation of the anatomy. A suitable method of calculating the NCC quickly is by the use of a Fast Fourier Transformation (FT), for the simple reason that a correlation 12 in the spatial domain corresponds to simple multiplication in the frequency domain.

After the relative translation has been calculated, the rotation of the patient on the treatment table can be determined. An FFT is particularly preferred at this stage as translation and rotation decoupie in frequency space. The rotation can be determined by calculating two new images which contain the absolute values of the fully transformed images from Section 5. One of these images is then rotated within a prescribed range and the NCC calculated by correlating the rotated image with the second image.. The maximum of the NCC then provides the rotation of the anatomy.

Figure 5 shows the Bandpass filtered image of Figure 4. The Bandpass filter responds to radiation field edges as well as to bony details and additional processing (section 6 above) is added to remove the borderline of the radiation field. (Figure 8). Bony ridges are expressed as negative values of p(xy). For an additional enhancement positive pixel values are set to zero and negative pixel values are inverted to positive values (Figure 8). The resulting image contains the position of bony ridges which is exactly the information necessary for image registration. Because of the elimination of the slowly varying grey-levels, the reduction of the air bubbles to their borderlines, and the removal of the field edge it is often possible to calculate the translation between the two images without choosing a R01 in the Reference Image, thereby using as much image information as possible.4 As a measure of similarity the Normalized Cross Correlation is used. For a fast implementation the calculation is done in frequency domain by applying a Fast Fourier Transformation. The correlation function cN(ab) is given by

FT -' (FT (p, (x, y)). FT (p, (x, y))) P, (X, Y) P, (X, Y) 1 2 xty Xfy where p'(x, y) and P2(Xl y) are two Bandpass filtered Portal images FT and FT1 mean forward and backward Fourier transform, respectively, means complex conjugation and the denominator is normalizing the correlation coefficient (CN is one for a perfect match of two identical images). The 13 maximum of the correlation function represents the best match of the Reference Image with the Portal image and delivers thereby the translation (a,b) of the patient on the treatment table. By calculating the complete correlation function it is possible to detect even large translations of the patient reliably which are most important for interrupting the treatment. Another virtue of the Fourier based method is that one is able to find the global maximum of the correlation function safely. There is no danger to stay at local maxima as it is possible using maximum searches in spatial domain like hill-climbing or Downhill simplex-method [231. By combining the translation of the patient (ab) with the translation of the radiation field (Ax,A) which has been calculated in 11.2. one gets finally for the translation (tx.ty) of the patient relative to the radiation field:

tx = a-A, ty = b-Ay If this translation exceeds a threshold (for example 5mm) the treatment is interrupted.

It will clearly be advisable to generate an alarm of some sort if the patient positioning error is above a certain preset threshold. For example, it could be an audible alarm or it could constitute an automatic deactivation of the radiation generator. Alternatively, a servocontrolled patient table could be provided such that after deactivation of the radiation field the patient can automatically be re-positioned according to the detected error. The eptire process is then preferably repeated with a new sample image.

Figures 9 and 10 show the correlation co-eff icients of the two greylevel images before and after processing to reveal the internal structure. Thus, Figure 9 shows'the correlation co-efficients of the bare grey-level images. It can be seen that whilst a clear maximum exists, the gradient to the maximum is relatively shallow meaning that the precise position of the maximum cannot. be determined accurately. Figure 10 on the other hand shows the correlation co-efficents of the two structure based images, and whilst there are many maxima in the image, one maxima stands out as being the clear best correlation. That maxima is sharply defined and it is a relatively simple matter to determine its position accurately.

14 RESULTS The above described algorithm has been tested with Portal images received from the Department of Radiation Oncology at the William Beaumont University Hospital in Royal Oak, Michigan, USA. 1139 pairs of Portal images of different anatomical regions (pelvis, lungs, head/neck region) were registered. The matching results are summarised in the following table:

Anatomical Region Pairs of Images Successful Matching Success Rate Pelvis AP/PA 544 536 99%.

Pelvis lateral 325 316 97% Lungs 144 137 95% Head nasal 51 51 100% Head temporal 39 39 100% icular region 36 33 912 1% r LAu 1139 1112 980 TO 0 tal The success of the matching was controlled visually by creating a difference image of the anatomical structure before and after the matching (Figure 12, Figure 13). The bony structures of the Reference Image are white and in the Portal image they are black. A human observer immediately realises when the matching does not perform correctly.

In addition, we examined the matching results when the ReferInce Image was either a DRIR or an image created interactively by drawing lines on top of bony ridges. We tested the matching of such a Reference Image with 17 Portal images from the pelvic region and were in every case successful. This was even the case when the image quality of the Portal image was especially poor or the Portal image contained artifacts like dosimeters. Figure 14 and Figure 15 show again two difference images before and after the matching of bony structures where the Reference Image is a highpass filtered DIRR.

An average success rate of approximately 98% is very encouraging and shows the superiority of the structure based matching algorithms to grey- level based methods.

There was one sequence of Portal images (lateral images from the pelvis) where the algorithm did not produce reliable results for the shift of the patient. The reason was that the relevant anatomical structures which had to be matched were lying very close to the edge of the radiation filed (the radiation field was rather small). By calculating the Bandpass filtered images and subsequent shrinking of the radiation field (as described above in 11.4), the relevant anatomical information was either completely or at least partially removed. As a consequence, even interactively it was difficult to identify the shift of the patient reliably. Therefore, this image sequence was not classed as a failure for the on-line algorithm. It simply emphasises the fact that the user has to examine the Reference Image carefully before the actual treatment starts and make sure that the relevant anatomical structures wich should be matched are inside the radiation field. An extreme example would be Portal images where no anatomical structures are present and where the above described algorithm would clearly be unable to calculate meaningful results. In addition, the user has sometimes to choose carefully an appropriate Region of Interest in the Reference Image to suppress disturbing features which should not be matched.

For Portal images of 256 x 256 pixels the algorithm takes approximately 1. 8 seconds on a SPARC-20 workstation which shows that the algorithm can indeed be used for a fast on-line control in the clinical environment.

Modifications to the above process are clearly possible. Forexample, median filters for the removal of point errors are not necessary as other methods can remove point errors reliably. Equally, segmentation of the radiation field area can be performed via other mathematical methods. Indeed, different methods of matching the radiation field areas can be used whilst still employing the second aspect of the invention as set out initially in this application. Noise suppression is not necessarily carried out by a Gaussian filter, and the block filter is not essential for removal of slowly varying grey-levels. Other high pass filters could conceivably be used. Likewise, the edge of the radiation field can be removed via a variety of methods.

13 16 The precise mathematical method for correlating the two structurebased images is equally immaterial to the present invention. The calculation can be done in the spatial domain, in the frequency domain or via a multiscale method, as can the rotation. The similarity function forcalculation of the correlation can be exchanged to alternatives to an NCC. What is important is that the translation and the rotation are determined on the basis of the structure based images, in order to convert the poor correlation of Figure 6 to the straightforward correlation of Figure 7.

The order of the above described image processing steps is not mandatory. For example, noise suppression can be combined with determination of the translation in Fourier space.

All such variations are naturally intended to be included within the scope of the present invention.

11 17

Claims (15)

1. A method of analysing a radiographic image comprising the following steps:

   (a) preparing a reference image; (b) preparing a sample image; (c) filtering at least the sample image to remove low frequency variations; (d) selecting a region of interest (R01) within the reference and sample images; (e) correlating the ROI of the reference image with the ROI of the sample image; (f) determining the relative displacements of the reference image and the sample image.
2. 2. A method according to claim 1 wherein the region of interest (R01) is be a selected field of the image.
3. 3. A method according to claim 1 or claim. 2 in which the sample image is an image obtained by:

   (a) preparing a reference image; (b) acquiring a portal image; (c) preparing a reference radiation field image and a sample radiation field image from the reference image and the portal image respectively by, in order, (i) defining a threshold value and (ii) setting pixels within the respective images to a dark or a bright state if darker or brighter than the threshold value; (d) optionally, comparing the reference radiation field image and the sample radiation field image to determine at least one of their relative rotation, translation and scaling; (e) subtracting one of the radiation field images from the other thereby to obtain a difference image.

   18
4. 4. A method according to claim 3 in which the sample image is shrunken prior to processing.
5. 5. A method according to claim 4 wherein the shrinking step is carried out by setting each pixel to a dark state unless all four neighbours are bright.
6. 6. A method according to claim 4 in which the sample image is shrunken by more than one pixel.
7. 7. A method according to claim 6 in which the shrinking step is c arried out by setting each pixel to a dark state unless all four neighbours"are bright, the neighbours being those pixels separated from the pixel concerned by a preset distance.
8. 8. A method according to claim 7 wherein the distance is between 15 and 30 pixels.
9. 9. A method according to any preceding claim wherein the reference image is one of a Portal image from the first treatment, a pretreatment test image, or a digitally reconstructed radiograph.
10. 10. A method according to claim 1 or claim 2 wherein the sample image is a Portal image taken during a radiotherapy treatment.
11. 11. A method according to any preceding claim wherein the sample image is captured electronically and analysed during the therapeutic radiotherapy treatment.
12. 12. A method according to claim 11 wherein data is captured during the first 10 seconds.
13. 13. A method according to claim 11 wherein data is captured during the first 5 seconds.
14. 14. A method according to any one of claims 11 to 13 wherein the radiotherapy session is automatically ended if inspection reveals a 19 significant error.
15. 15. A radiotherapy treatment apparatus comprising means for digitally capturing a Portal image and means for analysis of that image according to any preceding claim.

    1 11 1 1