PROTON RADIOGRAPHY WITH THE PIXEL DETECTOR TIMEPIX

This article presents the processing of radiographic data acquired using the positionsensitive hybrid semiconductor pixel detector Timepix. Measurements have been made on thin samples at the medical ion-synchrotron HIT [1] in Heidelberg (Germany) with a 221 MeV proton beam. The image contrast is generated by the energy of charged particles imparted after passing through the sample. Experimental data from the detector were processed for the establishment of the energy loss of each proton using calibration matrices. The interaction point of the proton on the detector was determined with subpixel resolution by a model fitting of the individual signals in the pixel matrix. Three methods have been used for the calculation of these coordinates: Hough transformation, 2D Gaussian fitting and estimate of the 2D mean. Parameters of the calculation accuracy and the calculation time are compared for each method. The final image was created by the method with the best parameters.


Introduction
Radiography is often associated with the X-ray radiography, where the image is formed by an attenuation of the beam after passing through the sample.The image is recorded on a radiologic film, a foil or a memory area detector [2].The process is similar to conventional photography, where negative is created according to radiation intensity.In contrast to this principle proton radiography can get the imaging information about one pixel with only one particle that has passed through this point [3].Therefore it is necessary to know the energy of the particle.Thus the principle of proton radiography is to determine the energy of each particle after passing through the analyzed sample.Hence proton radiography is suitable for imaging thin samples.The advantages of proton radiography are high sensitivity, high spatial resolution imaging of thin samples and the imaged objects receiving a lower dose [4].Although the absorbed dose by a single proton is significantly higher than that of a photon, the total radiation dose is reduced due to the need for only a single proton to display a single pixel [5,6].
The aim is to create a system of automatic image processing of proton radiography.Proton radiography, unlike traditional X-ray radiography, offers new possibilities and specific imaging thin samples, which can achieve greater sensitivity and spatial resolution.

Hybrid semiconductor detector Timepix
The Timepix detector [7] is a position-sensitive semiconductor pixel detector of the family of Medipix detectors developed at CERN [8].Timepix, based on the previous device Medipix2 [9], Timepix provides perpixel energy sensitivity.The detector consists of a radiation sensitive semiconductor sensor bump-bonded to the Timepix ASIC readout chip.The sensor can be of different material (Si, CdTe, GaAs) and thickness (e.g.300, 700, 1000, 1500 µm).The detector provides an array of 256 × 256 pixels for a total of 65 536 independent channels.Size of one pixel is 55 µm and the full sensor size is 14 × 14 mm (1.98 cm 2 ).Each pixel is connected with its analog and digital signal chain including amplifier, discriminator, counter and digital integrator [7].

Timepix principle
Ionizing particles produce in the sensor a cloud of charge which is collected as in a standard semiconductor diode.The collection of charges undergoes charge sharing (due to diffusion) and results in a signal containing several pixels forming a so-called cluster of pixels [10].The distribution of charge around the impact point has the character of a Gaussian.Data are recorded in frames which are saved to individual file frames which are time stamped.All data of events are saved to one frame during the selected time.Events which have been recorded at next selected time are saved to next frame.
The detector pixels can be independently operated in three ways depending on the mode in which the detector works [7]: • Medipix mode: Each event is counted in this mode.It can be used to count the event rate.
• Timepix mode: Value which is recorded in this mode corresponds to the time of interaction during the frame.
• Time over threshold (ToT) mode: The time registers the amount of charge which is given by the energy deposited.
Data which are measured by mode ToT is necessary to calibrate into energy.Calibration matrices are used for this calculation [11].Four calibration matrices are used A, B, C and T. The calibration coefficients are obtained by calibration procedure using discrete-energy X-rays from laboratory X-ray sources on XRF [11][12][13].The energy was calculated according to equation where X is the matrix with pixel values in energy in [keV], Y is the matrix which contains pixel values recorded in [ToT] counts and matrices A, B, C and T contain calibration values.

Description of experimental setup and radiographic data
Radiographic data used in this work were measured at the HIT ion synchrotron [1] in Heidelberg (Germany).
The sample was irradiated by almost monoenergetic proton beam with energy 221 MeV.A sample consisting thin foils which were stacked one over the other and a fly's wing were used a sample.The sample was irradiated in air.The Timepix detector which worked at ToT mode was placed behind the sample.Protons were registered by the detector after passing through the sample.Data from the detector were structured into 34 files.Each of these files contained around 1000 frames.The duration of one frame was 0.2 seconds.The value of ToT and coordinate of each pixel where the threshold had been exceeded were stored.The average number of protons events per frame was about 10 events per frame.The impact of proton on detector creates a round track with a diameter around 7 pixels as shown in Fig. 1.This event has as character so called cluster (see Fig. 1 and Fig. 2) that can be fitted by a 2D-Gaussian.The highest ToT values are located in the center of the track (Fig. 2).The position of impact on the detector can be accurately determined based on analysis and fit of all pixels of the event recorded.Sub-pixel resolution can be achieved as well.The total energy deposited to the detector by one proton corresponds to the sum of the energy values of all pixels in the cluster.

Data processing methods and their comparasion
Because all data that were contained in the files were not suitable for processing, it was necessary to select a suitable event record, to avoid unadvisable events and prepare the data for further processing.2D filtering was used for eliminating unwanted and false events.Nonzero pixels were calibrated using energy calibration matrices.Furthermore approximate coordinates of each protons impact have been found.An area of 11 × 11 pixels was cropped around these approximate coordinates.Because the diameter of the event record was about 7 pixels, it was probable that each pixel corresponding to this event record was included in this cropped area.Then the area was evaluated whether it would be suitable for further processing.If two or more event records overlapped at one cropped area, this area was not suitable for further work.It was necessary to accurately determine the location of protons impact for image reconstruction.Three methods were used for calculation of the precise location of the impact of the proton: Hough transform [14], 2D Gaussian fitting [15] and estimation of 2D mean [16].These methods were mutually compared and the reconstruction of image was performed from the results of the most suitable method.
Hough transform.The detection of the circle, described by equations (1), was used for assessing the center of the event by Hough transformation.The idea was in presumption that each event record had a circular nature.A matrix with data of the one frame was converted into a binary map.All contour points of event record were found by using the edge detector.
Coordinates of circles points which had the same radius as a radius of event records were calculated from the coordinates of each edge point.Eighty points of circle coordinates were calculated for each edge point.The value of one was added at each point, whose coordinates were calculated in this manner.Maximum of added values was subsequently calculated for each event.There is the impact point.Equations (1) describes the calculation of the circles coordinates with center at each edge points, where x 0 , y 0 are coordinates of edge points, θ is an angle from 0 to 2π.
Step of angle θ calculation was chosen up to 1/80 of the circle.Variables x(θ), y(θ) are coordinates of the new circle for one angle θ, r is the radius of the circle [14].
2D Gaussian fitting.Gaussian fitting was used for calculating center coordinates the impact of one Estimation of 2D mean.Estimate of 2D mean was used as the last method for determination center of events records.The energy values of events records were normalized, so that the total sum of all values was equal to one.First of all the normalized values were added up vertically.Then the results of these sums were multiplied by the respective horizontal coordinates.Finally these values were added up too.This result was one number only and it was the horizontal center coordinate.The same method was used for the calculation of the vertical center coordinates.

Comparision of methods
Methods have been mutually compared.Similar accuracy was reached using methods estimating 2D mean and 2D Gaussian fitting.Calculation of the estimate of 2D mean value was about 30 times faster than the 2D Gaussian fitting.Calculation using the Hough transform amounted to less accuracy, and it was also significantly slower than the calculation of 2D mean.
Comparing of the average Euclidean distance between the locations of particles impacts detected using each method and their standard deviations are shown by using each method in Table 1.It is evident from Table 1 that similarly accurate results were achieved using methods: 2D Gaussian fitting and the estimation of 2D mean.Mean Euclidean distance between of these methods was caused by mutual displacement in one direction.This shift was caused slight displacement of the resulting reconstructed image only.The representation of centers of events records using different methods is shown in Fig. 3. Histograms which are shown in Fig. 4 and 5 were created on the basis of the calculation Euclidean distances between the centers coordinates of events records detected using each method.

Image reconstruction
The final image was created from centers coordinates of each trace of impact whose calculation was described in the previous section.This image was created from the sum of the total energy.The value of the total energy which was given to the detector by a proton was calculated as the sum of all energy values of one event record.This value was assigned to the point whose coordinates have been detected as the place of protons impact.If one pixel of the final image was located for more events than one, the average value of these energetic records was assigned this point.Points without localized event are assigned zero value.These points appeared as holes in final image.Values for these points can be improved in the final image by extrapolating values from the surrounding area.
Just one proton can suffice for displaying a pixel.Each proton passed through the sample carries the information about a point of the sample.The final image can be displayed magnified due to the subpixel resolution of located impacts coordinates.The picture created from all the 33 000 frames was refined at double zoom.Localized have been around 340 000 events.This corresponds to an average of 1.3 protons per pixel.An example of such image is shown in Fig. 6.

Discussion and Conclusion
The resulting images, which were reconstructed from data obtained using the methods of 2D Gaussian fitting and estimation of 2D median values were essentially identical.To display a single image pixel was sufficient to locate the coordinates of only one proton impact.Because here we do not measure the decline of are irradiance detector, but directly the decrease of the energy particles is achieved by lower radiation load samples compared with other radiological imaging methods for example x-ray radiography.Experimentally it was found that it was ideal to enlarge final image twice.In this view, it may well calculate the value of the missing pixel, because in the vicinity of missing values there are many pixels designed so as to have one missing value extrapolated.The algorithm was created for a complete data processing.The data were processed automatically starting from the records of the detector and the resulting final image.Another advantage of proton radiography in addition to lower radiation burden is the availability of application accelerators.
For the future it is planned to apply radiography using charged particle accelerators Van de Graaff HV2500 -IEAP CTU [17] and cyclotron U-120M -NPI ASCR [18], which would also examine pokes energy charged particles after passing the samples.Also, the plan is to try and apply electron radiography to Microtron MT25 -NPI ASCR [19,20].At the electron radiography should be an option to display an image using electron scattering in the sample.Here, the dispersion is measured and determined by using the aperture, which should be placed on the Fourier plane between magnetic lenses respectively quadrupoles [21].It is also planned to develop a new automated data processing systems from new applications and other optimization of automated processing system previously available data.

Figure 1 .
Figure 1.One frame with cropped selection.The left figure shows one frame of radiographic data.This frame contains 9 suitable events and few unwanted values which are recorded by X-rays.Three suitable events and three unwanted values are shown on the right figure, where it is displayed the cropped area from left figure.

Figure 2 .
Figure 2. Registration of a proton event in Timepix operated in ToT (energy mode).The color and height of bars represents the energy of one pixel of the detector.

c 3 c 4 2 − y − c 5 c 6 2 ,
proton as another way.It was performed for only cut selection around traces of impacts.Data were fitted to 2D Gaussian with equation f (x, y) = c 1 + c 2 exp − x − where c 3 a c 5 are coordinates of the center of 2D Gaussian.It corresponds to the location of a proton impact.Coefficient c 1 corresponds to dislocation, c 2 is amplitude and c 4 and c 6 corresponds to spread of the blob.

Figure 3 .
Figure 3. Centroid fitting for one event recorded.The calculated centroid coordinates are shown in the middle of this figure.Data1 shows centered calculated using the 2D Gaussian fitting, data2 shows calculated center using the estimation of 2D mean and a data3 shows the detected center using Hough transform.

Figure 4 .
Figure 4. Histogram Euclidean distances between centroids computed by Gaussian fitting and estimation the mean 2D.

Figure 5 .
Figure 5. Histogram Euclidean distances between centroids computed by Gaussian fitting and Hough transform.

Figure 6 .
Figure 6.Final radiogram.Proton radiography of a composite sample consisting of an assembled foil array and a fly's wing.The final image was made from coordinate locations of proton impacts using estimating the 2D mean.One pixel represents 22.5 × 22.5 microns.There was located on the average 1.3 protons per pixel.

Table 1 .
Comparison of euclidean distance between the location of protons impacts and its standard deviation.