Using DosiVox to Reconstruct Radiation Transport through Complex Archaeological Environments

The DosiVox programme is used to reconstruct radiation transport in complex depositional environments. Using two archaeological case studies from ancient Egypt, the burial environments for a selection of ceramic vessels are reconstructed using the DosiVox programme, allowing the simulation of the emission and transport of radiation throughout these burial environments. From this simulation we can extract the external dose rate of the archaeological samples, a measurement necessary for determine a luminescence age. We describe in detail how DosiVox can be used to best advantage at sites with complex depositional histories and highlight that DosiVox is a valuable tool in luminescence dating. This work illustrates that DosiVox is, at present, unparalleled in reconstructing a more accurate and detailed external gamma dose rate which can significantly improve upon simplistic scaled geometric models.


Overview
Luminescence dating, including optically stimulated luminescence (OSL) dating, is a valuable tool to archaeologists and geologists and remains arguably the best technique for the chronometric dating of archaeological or geological materials comprised of quartz and feldspar. In this manner, in archaeology, luminescence dating is widely applicable to pottery, minerogenic building materials and ceramic objects or sediments from archaeological contexts. While dating of recently excavated material can be routine (as in situ sampling can be done by a trained luminescence specialist aware of precisely what samples and measurements are required to achieve a luminescence age), working with more complex material, such as that recovered from museum contexts, can be more difficult. While working with recently recovered material suitable for OSL sampling is of course ideal, it is not always possible and in many cases OSL dating is not permitted, owing to strict archaeological laws which coincide with a lack of accessible laboratory facilities.
Archaeological materials housed in museums-most notably pottery vessels and ceramic objects-were often not collected with OSL dating in mind and in some cases were excavated and placed in collections well before OSL dating was even developed. As such, external dose rate (Ḋ ext ) calculations can be incredibly difficult to reconstruct for museum material. Sometimes, for example if no external sediment from the original burial environment is forthcoming, it can be impossible to accurately reconstruct. However in some circumstances we can find museum material which is more favourable to OSL dating, for example, when original depositional sediment is attached to the vessels or is associated in some other way with the material. However, even then, the assumption of a uniform burial environment in archaeological luminescence dating is one that has often proved difficult for archaeologists to reconcile.
With the recent development of the DOSIVOX radiation transport modelling software [1], it has now become possible to reconstruct radiation transport throughout complex burial environments, where one must consider natural radiation from other archaeological materials, associated geological materials, and structural elements of the environment. Such modelling can vastly improve the accuracy ofḊ ext determinations for both museum specimens, as well as for complex burial environments recently excavated.
In this paper, we outline the general process for building a DOSIVOX model of an archaeological site or other complex depositional environments. This model, along with the DOSIVOX software, allows us to calculateḊ ext from first principles and certain basic assumptions, using all the information we have at our disposal (assuming direct measurement is not possible). Thus, for complex depositional environments where direct measurement is not possible, the model represents the best possible estimate ofḊ ext using tools available at the time of writing.
In Section 1.2 we describe our approach to the measurement ofḊ ext starting from first principles. In Section 2 we begin by describing the application of the model to two complex sites from Early Dynastic Egypt where direct measurement ofḊ ext is not possible. We then, continuing into Section 3, compare the results to other possible methodologies for estimatingḊ ext . In Section 4.1 we describe how to construct a DOSIVOX model based on known information about a site and in Section 4.2 we describe our approach to running DOSIVOX, including ensuring the quality of the simulations.

Theoretical Background
In our work, we start by observing that, owing to their very short attenuation depth, neither alpha nor beta radiation originating externally to the sample contribute to the equivalent dose. In other words, we assumeḊ ext is entirely composed of gamma radiation.
With that assumption, three key questions need to be answered when modellingḊ ext : 1. what was the quantity (flux) of gamma radiation arriving at the target site (the location of the sample) per unit of time; 2. what fraction of that radiation was absorbed; 3. what dose of laboratory-delivered beta radiation would be equivalent to the absorbed dose of gamma radiation?

Grain Absorption and Equivalent Beta Dose
We begin by tackling the second and third questions, which we will find, in our case, to be less interesting than the first. To estimate the fraction of incident radiation that is absorbed by a grain of quartz, it would be possible to model the grains themselves using DOSIVOX. However, given the extreme difference between the grain size and the gamma attenuation depth (roughly 1:2500) it is highly accurate to model the grains as thin targets, as is the usual practice [2]. This means that we assume that grains are irradiated evenly. Thus, we omit grain models and use DOSIVOX to calculate the average absorption in a uniform material of the right composition, as we will describe below.
In answer to the third question, we note first that the relationship between the luminescence caused by a dose of gamma radiation and the luminescence caused by the same dose of beta radiation remains impossible to determine from first principles. It would be possible to measure it experimentally by irradiating first with a gamma and then a beta source for the samples in question; however, it is in fact well established that the relationship between the two does not vary from one sample to another, and that it is accurate to assume that they are equal (Reference [3]: 47).
Thus in answer to the second and third questions, it is both standard practice and accurate to assume that the dose is evenly absorbed and the equivalent beta dose is equal to the gamma dose.

Incident Dose
Answering the first question will be the primary focus of this work. In general, the radiation incident on the grains can be directly measured using in situ gamma spectrometry, estimated using the infinite matrix assumption (IMA, see below), or, as in this work, calculated using DOSIVOX. We begin by briefly discussing gamma spectrometry, then examine in more detail the advantages and pitfalls of the IMA, and then lay out the new possibilities when using a first-principles method like DOSIVOX.

Gamma spectrometry
Gamma spectrometry is the direct measurement of the incident gamma flux at the site using a portable measuring device, that is, a gamma spectrometer (or a dosimeter, a gamma spectrometer may require an additional assumption about secular equilibrium). If it is possible to take a gamma spectrometer or dosimeter to the sample site and measure the incident dose directly, then that is undoubtedly the easiest method to pursue. However, there are reasons why this may not possible, for example: • A portable gamma spectrometer is unavailable or cannot be transported to the site.

•
The site has radically altered since antiquity; thus, the modern measurement of the gamma radiation intensity may not reflect what it was in the past. • The context of the sample has already been changed during the current excavation, for example as a result of the removal of stratigraphic layers above. • The artefact of interest was excavated in the past and its find spot either cannot be accessed or has been destroyed, for example, by backfilling.

Infinite Matrix Assumption Plus Geometric Model
Let us now turn to the case where direct in situ measurement of the incident gamma flux is not possible. A second method of determining the gamma flux begins by recovering samples of all materials within the gamma attenuation depth of the sample and determining the intensity of gamma radiation emitted by these materials. This may be through direct measurement using a lab-based gamma spectrometer, or by determining the radioisotopic concentrations within each sample using a form of mass spectrometery (e.g., ICP-MS). The gamma emissions of all commonly discovered radioisotopes are well-measured and readily available [4]. Thus, if we know the radioisotopic content of any material or material mix we know the amount of radiation emitted by that material (per unit mass). Regardless of which method is used, the amount of radiation emitted by each material can be determined. This is a good start, but what we want to know is the amount of radiation absorbed by the sample (per unit mass). Happily there is one simple case in which one is equal to the other: the case of a homogeneous material of infinite extent. In which case, by symmetry all points in the material must absorb the same radiative energy, and by the conservation of energy, this absorbed amount must be equal to the emitted amount, which we know (given the elemental breakdown). In the case where the material is not infinite, but is nonetheless considerably larger in extent than the attenuation depth of the gamma radiation (typically ∼30 cm in archaeological and geological contexts), we may treat the material as being infinite, and give the absorbed dose as being the emitted dose: this is known as the infinite matrix assumption (Ref. [5]: 56).
The IMA is an excellent assumption if the sample is indeed embedded in a uniform material, as shown in Figure 1, but what if instead it is located at the boundary between two or more materials? Provided the materials can be reasonably approximated as occupying solid angles surrounding the sample, all radiation arriving at the sample only passes through one material before reaching the sample and the infinite matrix assumption still applies, weighted by the solid angles (see Figure 2). It is possible to prove this in some circumstances and we do so in Appendix A.
The convenient arrangement of materials in Figure 2 is unlikely to apply in a real case. However, in the absence of any other methodology we may assume that we can still use this method, approximating the solid angles that apply to each material as the volume occupied by that material in a sphere whose radius is the gamma attenuation depth and whose centre is at the sample, divided by the total volume of the sphere. The more complex the actual arrangement of materials, the more specious this approximation, as our results in Section 2.4 demonstrate. In the case of contexts with a simple arrangement of materials, this geometric approximation will frequently be acceptable (e.g., Reference [6]; see extended discussion in Reference [5]: 289). In the case of a complex archaeological context like a tomb, something more is needed.

DOSIVOX
The large hadron collider (LHC) at CERN causes particles to collide with such high velocity that they break up into a shower of other particles, which are seen in the detector. Because the effects of these particles in the detector are indirect it is necessary to model the origin, transport and detection of these particles and then use statistical methods to determine the most likely sequence of events (e.g., Reference [7]).
One such model created and used at CERN is Geant-4, which simulates the creation, transport through matter and absorption of high-energy particles. Researchers at the Université Bordeaux Montaigne realised that this is exactly what was needed to model the dose rate that a luminescence sample experiences in its burial location, and thus the DOSIVOX software was developed [1].
As will be discussed further in Section 4.1, using DOSIVOX starts with making a model which includes the geospatial arrangement of all materials in the burial environment, as well as their composition. The composition of the materials includes their bulk composition, that is, the elements that make up the bulk of their mass (which affects the transport and absorption of the radiation) and the radioactive trace elements which emit the radiation (which are specified separately). We note that as well as specifying the composition and distribution of the materials, it is also possible to specify the grain-size of the sample; however, we omit this as we are only interested in the gamma radiation.
The dose absorbed by the sample location is the integral over time of the energy deposited by every gamma photon emitted from every radioisotope in the model. Since this integral is high-dimensional and expensive to calculate, DOSIVOX uses Monte-Carlo methods to approximate the integral, effectively allowing a given number of gamma particles to be generated at random from within the model and measuring the dose at the sample location from each particle. Thus, the greater the number of particles, the more accurate the calculation, as discussed in Section 4.2. For more details concerning the way DOSIVOX works the reader is referred to the original paper [1]. In this work, we wish to give a practical description of how DOSIVOX was applied, in the hopes that it may be useful to others wishing to use DOSIVOX in complex archaeological or geological contexts.

Other Gamma Sources
There are of course other sources of gamma radiation that might have affected the effective dose rate. Assuming that the sample was not exposed after excavation, the most notable of these is cosmic radiation. In our work here, we use the effective and well-established method of calculating the cosmic dose rate as a function of depth and global position of an archaeological/geological site that is provided by the Dose Rate and Age Calculator (DRAC) ( [8]; following Reference [9]).

Overview
To illustrate how DOSIVOX can be used to reconstruct radiation transport through complex burial environments, we look at two case studies from ancient Egypt. This research was initially carried out as part of Hood's doctoral thesis, which produced the first ever OSL dates obtained for ancient Egyptian ceramics. The results of this dating programme can be found in References [10,11]. Owing to a current law in Egypt which prevents the removal of archaeological material for scientific analysis, this research was carried out on material excavated at the turn of the 20th century which was stored in museums outside of the country. As a result the burial environment and thus DOSIVOX model, was established through analysis of the original excavation reports, and when possible, more recent archaeological examinations (as was possible in one case, see below).

Case Study 1: Bêt Khallaf
The first case study looks at twenty-four ceramics sampled from the Upper Egyptian site of Bêt Khallaf. The site was excavated by British archaeologist John Garstang in the early 1900s, and the intact ceramic material recovered from the site was removed from Egypt soon after excavation and dispersed across three museums: the Ashmolean Museum in Oxford, the Garstang Museum of Archaeology in Liverpool and the University of Pennsylvania Museum of Archaeology and Anthropology in Pennsylvania (Penn Museum). The assemblage is typologically characteristic of the late Early Dynastic Period of Egypt (c. 2650 B.C.E.).
The Bêt Khallaf assemblage is a perfect candidate for testing and comparing gamma dose rates obtained using DOSIVOX to those determined using the IMA. This is firstly because many of the vessel types found in the tombs at Bêt Khallaf are diagnostic (thus providing a perfect platform for using absolute dating to anchor the relative ceramic typology to a known calendrical system), and secondly because Garstang kept and published reports of high standard for that era. These reports were used to construct the DOSIVOX model. In addition to providing details of the tomb architecture (tomb size, depth etc.), Garstang also gives the location of the pottery deposit within each tomb, allowing an accurate find spot of the ceramics to be determined (Reference [12]: 8-16).

Case Study 2: Tomb of Djer, Abydos
The second case study was carried out using four vessels from the Tomb of Djer at Abydos. Excavated by Flinders Petrie in the early 1900s, this Abydene tomb was the burial place of King Djer, the 3rd king of the first Egyptian dynasty (c. 3000 B.C.E.). After excavation, the ceramic material was transported to the Ashmolean Museum, Oxford. Similarly to the Bêt Khallaf material, the Djer material was well-documented for its time and Petrie discusses the finding of the vessels and their precise location within the tomb. The DOSIVOX model was constructed based upon Petrie's excavation reports but was also supplemented by newly documented evidence from the tomb, which has recently been reexamined by the German Archaeological Institute. The combination of both sets of information allowed us to reconstruct the deposition of the vessels with a high degree of certainty and renders the Tomb of Djer a prime candidate for examining gamma dose rates and comparing the results obtained using DOSIVOX and the IMA.

Comparing Gamma Dose Rates
In both ancient Egyptian case studies,Ḋ ext was determined using the U, Th and K concentrations from original depositional material found adhering to the vessels (under the rim and on a vessel handle). DOSIVOX was used to model the gamma dose rate received by the vessels from the surrounding burial environment (where the alpha and beta dose rates were determined using DRAC). The cosmic dose rate (Ḋ cos ) was determined by using information reconstructed using the original excavation reports and Google Earth (where one can easily identify the tombs from which the vessels came from and thus their elevation, etc). Using the archival material discussed above, it was possible to reconstruct the burial environment for the examined ceramics. In both cases DOSIVOX modelling was considered necessary as the depositional environment for both sites was considered complex, without straightforward sample geometry. Complex sample geometry in luminescence dating occurs when a sample is not obtained from a (more or less) homogeneous depositional environment. In contrast, a simple sample geometry is what we see in Figure 1.
In the case of Bêt Khallaf, the ceramics were not fully buried, but rather 'strewn' about in a presumably haphazard manner. The precise manner of their deposition is unknown (i.e., they could potentially be sitting upright, lying on their side, or halfway out of the tomb fill; as described in Section 4.2.1 these variables were tested to see their effect upon the determined gamma dose rate and have been presented below), but from the tomb reports we were confident that they were not fully buried and thus would be affected by complex sample geometry, consisting of the pot itself, the tomb fill/debris on the ground as well as the tomb walls.
In the case of the Tomb of Djer, Petrie recorded the exact find spot of these vessels in the tomb, and the unpublished work of the German Archaeological Institute further informed us that these vessels were positioned upright and buried up to their mid-point in the surrounding sediment of a chamber in the north west corner of the mud-brick tomb and covered by a mud-brick staircase. These vessels also contained contents which were a mixture of organic and inorganic material. All of these known details helped to reconstruct the complex sample geometry for these vessels and to make the resulting DOSIVOX model.
Through reconstruction of the burial environment our resulting DOSIVOX models were able to present a more realistic gamma dose rate measurement for the OSL sample than if we simply assumed 4π geometry or determined a 'best guess' geometry based on assumed percentage contributions from each material measured.
Tables 1 and 2 present the calculated gamma dose rates for the Bêt Khallaf ceramics (Table 1) and those from the Tomb of Djer (Table 2). For each set of ceramics the gamma dose rate for each sample has been calculated in one of three ways: firstly by using a 4π geometric assumption using the IMA; secondly using a scaled IMA based on a 'best guess' of the complex sample geometry and; thirdly by using the gamma dose rate modelled using DOSIVOX. In the case of Bêt Khallaf, a sample 'best guess' of 2π is displayed. This was derived by following the principles outlined in Section 1.2.2, as the samples were known to be on the floor of the tomb (i.e., fill below and air above); in the case of the Tomb of Djer, the pots contained organic residue and were more than half buried in fill. Thus, a 'best guess' of 3π fill and 1π pot contents is displayed.

Discussion
We can see from the results presented in Tables 1 and 2 that the gamma dose rates obtained by using DOSIVOX differ from those obtained by assuming 'best guess' scaled geometries. In the case of Bêt Khallaf the outcome saw DOSIVOX produce a gamma dose rate far more in keeping with an assumed 4π geometry than in our 'best guess' scaled geometry, illustrating how DOSIVOX can provide a result which is not subject to inbuilt biases and assumed geometries determined by the luminescence practitioner.
It can, of course, well be argued that the 'best guess' was in fact 4π, because the roof of the tomb was of similar composition to the floor and the dry air of the tomb would only slightly attenuate the radiation. We do not dispute this; we merely observe that all of these sweeping arguments and assumptions can be replaced by a more rigorous procedure using DOSIVOX. One benefit of this rigorous procedure is that there are many hidden assumptions in the 4π IMA estimate that DOSIVOX makes explicit and testable (position of the pot, etc). More importantly, however, it is only having run DOSIVOX (i.e., with hindsight) that we can know that the 4π IMA assumption was acceptable. Even though it might seem a reasonable assumption about the transport of gamma radiation a priori, with DOSIVOX the transport and absorption of gamma radiation is modelled in a way (using Geant-4) that has been extensively tested against highly accurate particle physics experiments; thus we can be much more confident in its results and we can eliminate at least one big assumption from our chain of reasoning.
In the case of the Tomb of Djer the DOSIVOX gamma dose rate was significantly different to the IMA and 'best guess' geometry dose rates. The complex depositional environment and the array of different elements affecting the dose rate was such that the gamma dose rate achieved using DOSIVOX modelling varied significantly enough from assumed geometries that final dose rate determinations varied considerably (which in turn affects the final luminescence age calculations). (It is beyond the scope of this paper to present the final age calculations for this material. If interested, the reader is referred to References [10,11], noting that final ages therein obtained using DOSIVOX were in good agreement with associated archaeological and historical chronologies.) In simple geometric cases, for example, when a sample is known to have full, 4π geometry, DOSIVOX may be too significant an undertaking to apply to every luminescence sample and indeed will, in all likelihood, give a similar result. However, the more informative analysis and more accurate gamma dose rates achievable by DOSIVOX modelling in more complex burial environments, for example, that seen in the Tomb of Djer, illustrates how using the programme can improve luminescence age determination. In the case of complex sample geometries, it is therefore recommended that DOSIVOX modelling be incorporated into luminescence work to improve dose rate calculations and, in turn, age determinations.
One limitation of this research is that from the tomb reports we know that the ceramics at both sites were found as caches, that is multiple vessels were found in the same deposit. However, it was impossible at either site to determine the order or precise location of the ceramics relative to one another. In an ideal situation, if this information were recorded at the time of excavation, it would be possible to improve upon the DOSIVOX model even further by adding in each individual ceramic into the model so that the gamma dose rate determined for each vessel would also factor in the gamma dose being emitted by the surrounding vessels. This could be easily achievable for those wanting to use DOSIVOX modelling of complex deposits found in current excavations or fieldwork situations.
DOSIVOX is a valuable new tool when dealing with samples that have complex geometry and where in situ gamma measurements are not able to be made. Using DOSIVOX in conjunction with current projects could yield gamma dose rates far more accurate than those achievable through a simple incorporation of the IMA. This is because, in most cases, all the elements required to construct a DOSIVOX model would be readily available and, as a result of the more detailed model and more limited (and explicit) assumptions, a higher degree of confidence could be placed upon such a model, thus potentially yielding more accurate gamma dose rates than the IMA alone. We recommend that DOSIVOX modelling is routinely carried out at locations where in situ measurements are not available.
DOSIVOX modelling can be used in both helping to determine the gamma dose rate at recently excavated sites, those sites excavated in the past (i.e., where the resulting material has been stored in museums) and for both geological and archaeological situations. We are of the opinion that the incorporation of DOSIVOX into luminescence research can be of significant value in a number of situations, and believe that more accurate gamma dose rates can be achieved using DOSIVOX in cases of complex sample geometry. As such, DOSIVOX has the ability to greatly improve the reconstruction of dating profiles at geological and archaeological sites worldwide and assist in more accurate chronologies being achieved.

Constructing a DOSIVOX Model
A DOSIVOX model broadly consists of two parts: 1. a model of the composition and layout of the materials in the burial environment, and 2. a similar but more detailed model of the target site (the place where the sample is drawn from).
Before discussing how those parts are created it is helpful define the following terms: • Material-the model is divided into regions that are filled with materials that are broadly uniform in composition: for example, a layer of sand, or a mud-brick wall. • Component-materials are described by their density, their water content and the relative quantities of their chemical components (e.g., SiO 2 or TiO). • Voxel-the whole burial environment is broken up into a grid of cubes which are called voxels (effectively "3D pixels"). Each voxel contains one material. • Subvoxel-in the mode of DOSIVOX used in this paper, one voxel in the model may be broken down into a fine grid made up of subvoxels. This voxel is referred to as the 'sub-voxelized voxel'. The sub-voxelized voxel is used when one part of the model needs to be described in more detail.
In this work, we have the case of a ceramic vessel in a tomb. The tomb itself is described at lower resolution using the voxels; one of these voxels (the sub-voxelized voxel) contains the vessel, which is described at higher resolution using subvoxels. This is elaborated below in Section 4.1.1.

The Model Layout
The volume that is being modelled (e.g., a tomb) must first be reconstructed, and then represented in terms of a 3D regular grid of voxels. Reconstructing the volume may be done using anything from field reports and excavation maps, to the latest 3D mapping tools like GIS. Here we discuss our own work using excavation reports from the early 20th century.
As discussed above, a major benefit of the Bêt Khallaf and Djer material is that Garstang's and Petrie's excavation records [12,13] were of excellent quality for the early 20th century and they allow us to reconstruct with a reasonable degree of confidence the exact burial environment of each vessel, which is further strengthened by finding parallels for similar temporal contexts and recent excavations (i.e., in the case of the Tomb of Djer).
Regardless of the source of information about the burial environment, what must eventually be constructed is a series of evenly spaced layers (i.e., horizontal sections) each of which is composed of an identical regular grid of voxels. In our model, the voxels were sized 30 cm × 30cm × 15 cm, with the total number of voxels varying with the size of the tomb being modelled. Thus, once the tomb shapes and sizes were reconstructed, each tomb layout was replicated using voxels in a three-dimensional grid.
Every voxel interior in the grid must be labelled with a number identifying a particular material. One of the voxels on one of the layers must be chosen as the sub-voxelized voxel-in our case, the location of the vessel of interest (as far as could be known, as we discuss below). Within the sub-voxelized voxel the whole process is repeated again on a finer scale. The voxel is divided into layers, which are divided into rectangular cells. Each cell, that is, each subvoxel, is labelled by a material representing a number (this process allows a higher spatial resolution in a given section of the tomb (i.e., the section that contains the vessel, without the cost of having such a high spatial resolution across the whole tomb ( [14])).
A variety of possibilities are available when constructing the sub-voxelized voxel. As discussed in Reference [14], there are many 3D scanning tools which may be used to make a 3D image of the vessel, which can be converted into a series of subvoxel layers by, for example, associating colours with materials. In our work, we constructed the sub-voxelized voxel manually from drawings of the vessels. While we were careful to preserve the overall scale and thickness of the vessel, we determined it was unnecessary to precisely specify the shape of the vessel, owing to the long attenuation depth of the gamma rays. Accordingly, each jar-shaped vessel was simulated using a simplified model of the vessel composed of 1 cm × 1 cm × 1 cm sub-voxels forming a regular rectangular shape (and for the bowls and the pot stand other approximate shapes were used). Figure 3 presents a set of schematic representations of the DOSIVOX simulation model used for vessel X5472 in tomb K2 at Bêt Khallaf. Figure 4 presents a similar set for vessel X5476 in the tomb of Djer. Supplementary Materials Video S1 presents a representative 3D visualisation of a ceramic vessel in a tomb.

Describing Materials
Having labelled all voxels and subvoxels as containing a particular material, we now need to define these materials. For example, the materials of the tombs at Bêt Khallaf included air (within the tomb chamber), a sand-based fill (on the floor of the tomb), and gebel-Egyptian bedrock-which comprised the walls and roof of the tomb.
To define a material, we require its density, its bulk chemical composition, the concentrations of radioisotopes, and the water content. As with the tomb dimensions, the methodology for determining this information for each material will vary depending on the project. In an ideal situation, a large sample would be obtained of each material, which would then be sent away for analysis using, for example, a form of mass spectrometry. However, in most cases, all that is possible is to incorporate all available information, and to make assumptions where necessary. The assumptions themselves can then be tested, as we demonstrate below.
Below we present a complete list of all the materials used in this project, the methods used to analyse them and necessary assumptions. • Air-assumed to be typical atmospheric concentration, with water content based on a specified temperature. • Fill at Bêt Khallaf -bulk chemical composition and radioisotopes determined from analysis of sediment attached to the outside of the vessels; density assumed to be that of loose sand.
• Gebel at Bêt Khallaf -chemical composition and water content assumed to be the same as the fill; density assumed to be that of sandstone; • Fill at Djer-radioisotopes determined from analysis of sediment attached to the outside of one vessel; bulk concentrations assumed to be similar to those at Bêt Khallaf; • Clay (vessel fabric)-measurement of the bulk composition is discussed below in Appendix B.1; radioisotope concentrations were measured separately for every vessel using ICP-MS analysis of samples extracted using the minimum extraction technique [15] (the sub-voxels within the vessel were classified as two separate materials-'Clay' and 'ClayBase'-so that the average dose rate was recorded for the base (i.e., the OSL sample location) and the rest of the vessel separately. However, these two clay sub-voxel types were identical in their composition). • Pot contents at Bêt Khallaf -vessels at Bêt Khallaf were assumed to be half-filled with fill. • Pot contents in the Tomb of Djer-the vessels in the Tomb of Djer contained the organic/minerogenic residue of their original contents; this residue had been analysed as part of another project [11,16,17]. • Mud-brick walls in the Tomb of Djer-no direct samples or measurements were available of the mud-brick composition. Accordingly, the composition of the mud-brick was assumed to be similar to that of Nile silt clay.
The full chemical breakdowns of each material that resulted from these analyses and assumptions are detailed in Appendix B.

Specifying Radioisotope Concentrations
While the bulk chemical composition is specified only once per material, the radioisotope concentration is specified on a per-voxel (or per-subvoxel) basis. This allows the same material to have different radioisotope concentrations in different places. Accordingly, as well as a grid specifying which material is in which voxel (and another for subvoxels), we must also provide grids specifying the radioisotope concentration in each voxel/subvoxel. DOSIVOX calculates the dose from each radioactive element separately, and so a separate grid must be provided for each one.

Running and Resolving DOSIVOX
Once all the necessary data is assembled, it is time to input that data into DOSIVOX and then use it to calculate the dose rates. Running DOSIVOX consists of the following steps. For the first step readers are referred to the DOSIVOX documentation [14]. For carrying out the second step, DOSIVOX provides a graphical tool. However it is also perfectly possible to manually construct a pilot file based on the included files. In this work, a template pilot file was manually constructed for each sample studied. From that point, an automated tool called CODERUNNER [18,19] was used to 1. set the radioisotope concentrations in the pilot file (for U, Th and K in turn); 2. run steps 3 and 4 above repeatedly, increasing the number of particles until the dose rates converged (that is, until the dose rates are known to the required precision); 3. parallelize the process by running multiple identical simulations concurrently, as detailed towards the end of this section.
Because of the Monte-Carlo methods used by DOSIVOX, there is no concept of time in the simulation. What is available in the DOSIVOX output data is the total dose emitted (E γ ) and the total dose absorbed (D γ ) by each voxel and subvoxel during the course of the simulation. In order to convert that total absorbed dose into a dose rate (Ḋ γ ), we must scale it using the emitted dose rateĖ γ (a known quantity given the concentration of the particular radioisotope, obtainable from, for example, Reference [4]). Thus, the absorbed dose rate from a given voxel or subvoxel iṡ This conversion was carried out automatically by the CODERUNNER tool. The output D γ is provided voxel by voxel and subvoxel by subvoxel. However, it is also provided as an average per material. As mentioned earlier in Section 4.1.2, we define the material of the base of the vessel (where the samples were taken in every case) separately to the walls. Thus in every case we are calculating the average dose rate in the base of the vessel.
Given the complexity of the model and the relatively small volume occupied by the target site (the base of the vessel), a very large number of particles were required to achieve a sufficiently low uncertainty in the dose rates. Since the simulations were carried out on a laptop (albeit one with an up-to-date hyperthreaded 4-core CPU) total CPU time was limited. Nonetheless, even with the limited simulation time available, almost 550 simulations were carried out for this project (using approximately 5000 CPU hours). High-resolution simulations used 8 million particles; lower-resolution simulations used 2 or 4 million particles. For each case, six individual simulations were run in parallel (i.e., the same case, and hence the same input file, was run in parallel). As DOSIVOX uses a Monte Carlo method it was possible to take the mean of the six simulations and combine these to increase the simulation resolution to 48 million particles per high-resolution case, that is, six parallel simulations × 8 million particles (Pers. Comm. L. Martin, 2015).

Testing Assumptions in the Reconstruction
As touched on above, the information provided about the burial environments were incomplete, and several assumptions had to be made when constructing the DOSIVOX models. A DOSIVOX model allows us to obtain the best possible estimates given available sources of information. However, importantly, it allows us to assess the size of the error that may result from these assumptions. Thus, DOSIVOX allows us to be quantitative about the effect of missing information from the archaeological or geological record.
In this section we detail several of the necessary assumptions that were made for this project, and assess their impact on the accuracy of our results. These are assumptions concerned with: 1. The orientation of the vessel: the pot can be either standing (upright or upside down), or lying on its side. Since the OSL sample is drilled from the base in all instances in this project, the orientation of the pot in its burial environment will affectḊ γ (i.e., was the sample location submerged in fill or exposed to air?) 2. Whether standing or lying on its side, the pot can be buried from 0% to 100% in fill, with the remainder surrounded by air (or other material). 3. The vessels can be filled from 0% to 100% with either fill or other contents (e.g., organic bulk residue). 4. The density of the vessel's clay matrix. 5. The density of the fill. 6. The vessel's proximity to the tomb walls, that is, its location within the tomb. Was it in the centre, corner, or along a side of the tomb? 7. The moisture content of the burial environment.
With regard to the position and burial situation of the vessels found within chamber tombs (i.e., where the tomb is not entirely filled up with fill, such as in the cases of Bêt Khallaf and the Tomb of Djer), parallels at recently excavated sites can help inform of the most likely way in which the vessels were found buried. Figure 5A provides a schematic representation of how a ceramic vessel cache may look in an Early Dynastic chamber tomb (based on Hood's personal observations in the field). As can be seen, the vessels in such caches can be positioned in a number of ways and Figure 5B-E show the four most likely of these. Owing to time constraints caused by the significant simulation time needed to test each case, this project limited itself to testing the first five assumptions listed above at high resolution, as these were deemed the most likely to affectḊ ext . Of the final two assumptions (assumptions 6-7), assumption 6 was considered, but at low resolution (the attenuation of radiation through air with low moisture content is such that it should not matter where the vessel is located as it will, in theory, still receive radiation from all four walls of the tomb). Assumption 7 was not tested as the literature consulted seemed unanimous that the correct approximate value for water content within the burial environment was 3% ± 2% ( [20][21][22][23]). One situation that can not be tested for in these case studies is the proximity of other vessels to the case vessel. It was impossible, for lack of information, to determine the order or arrangement in which the vessels were located in the tombs. Thus, each vessel was modelled to be the only vessel in the modelled tomb (see Section 3).
A full description of the tests carried out for these seven assumptions and the effect of their changes onḊ γ is given in Table 3 and illustrated in Figure 6. Figure 6. Graph depicting the effects of changes to the assumptions made in the DOSIVOX modelling. The change in the final age in years as a result of changing each assumption is plotted (using internal dose rate and equivalent dose rate data from Reference [10]). For each assumption we change, we display the effect on the age as determined from the coarse grain fraction (CG) and the fine grain fraction (FG) from the sample (NB a final age would result from combining these two values). The inner set of larger bars represent the magnitude of the error in years (at the 1σ level) resulting from the uncertainty in the DOSIVOX result for the reference case. It can be seen that changing all assumptions results only in errors which are comparable to or less than the typical error in the DOSIVOX simulations (which is itself small compared to the overall error). This means that there is insufficient evidence to reject the null hypothesis (the null hypothesis being that there is no effect of changes, or more precisely, that the effect of changes is smaller in magnitude than the uncertainty resulting from the DOSIVOX simulations). The outer set of larger bars represents the magnitude of overall uncertainty of the measured sample age in years. We can see that the DOSIVOX-derived uncertainty is itself only a minor contribution to the overall uncertainty, and that since the null hypothesis states that the effect of changes in assumptions are small, i.e., undetectable, with respect to the DOSIVOX-derived uncertainty, we can therefore conclude that the effect of changing the assumptions is negligible.
Overall it was found that changing these assumptions only had a weak effect (i.e., if there was an effect its magnitude was not statistically significant). This somewhat unexpected result is certainly worthy of further study. It is likely to be a consequence of the particular circumstances being considered here (specific dosimetries and so on). However, it may be hypothesised that it could be due to the characteristics of the emission spectra of the radionuclides. The emission spectra of U and Th are dominated by lower-energy, short-range gamma rays. Since the concentrations of U and Th are larger in the clay than in the surrounding fill, it might be expected that the dose rate from these radionuclides is dominated by the dose rate coming from the clay of the vessel (effectively the internal gamma dose rate) and will thus be unaffected by most of the assumptions tested here. By contrast the K spectrum is dominated by one high-energy, long-range gamma emission, for which the attenuation rate in any of the components of the model of the tomb will be low. Therefore, whatever the orientation of the vessel or its immediate surroundings, the absorbed K dose rate will effectively be an average of the tomb as a whole, i.e., it will be affected by whatever material is within 30 cm of it, or 30 cm of the material which is separated from the pot by air (as the gamma radiation will travel directly through air with low moisture content). Thus, even though the incomplete archaeological literature produced unavoidable uncertainty in the model, this did not give rise to significant uncertainty in the calculatedḊ γ . Table 3. The assumptions tested by DOSIVOX, illustrating the possible scenarios/assumptions that could be made (left column) and the outcome of the DOSIVOX modelling for each assumption (right column).

Modelling Unknown
Effect on Modelled Dose Rate of the Bêt Khallaf Assemblage 1. The pot can be (1) standing (upright or upside down), or (2) lying on its side. Since the OSL sample is drilled from the base in all instances, it also matters whether the sample location was on the upper or lower side of the deposited vessel (i.e., was the sample location submerged in fill or exposed to air)?
Two scenarios were considered: the pot standing on its base or lying on its side. These were considered the two limiting cases and as no measurable change was observed, it was concluded that the orientation of the pot and the location of the sample had no significant effect. To further examine the seeming minimal effect of vessel orientation on the measureḋ D ext , it would be possible in the future to determinė D ext for a set of samples taken from different locations on the same vessel.
2. Whether standing or lying on its side, the pot can be buried up to 100% in fill, or be unburied and thus exposed to air up to 100%.
Two scenarios were considered: the vessel lying on the surface of the fill or buried ∼15 cm, which for bowls means they are completely buried and for storage jars, they are buried up to half their height. Again, no measurable change was observed.
3. The vessels can be filled from 0% to 100% with either fill or other contents (e.g., organic bulk residue).
Two scenarios were considered: the vessel was filled with fill or it had no fill inside it. Again, no measurable change was observed.
4. How does the density of the vessel's clay matrix affect the externalḊ?
Although we were not able to measure the clay density of all vessels, a single sherd (X5460) was measured as having a density of 1.4 g/cm 3 . In the visual analysis of the clay fabric that accompanied each of the vessels the density of the vessels was determined as being 'high', 'medium', or 'low'. Independently of the density measurement carried out upon X5460, the vessel this sherd came from was deemed to have a 'low' density. Thus, based on this information, we were able to extrapolate that the density for vessels deemed to have a medium density was likely to be ∼1.6 g/cm 3 , and for vessels with a 'high' density, ∼1.8 g/cm 3 . Two scenarios were considered-the clay density was either set to its known value or to a reference value of 1.8 g/cm 3 . No measurable change was observed between the two scenarios and thus it was concluded that although the clay densities are not known precisely for each vessel, variation within the likely range of values will not affect the result. Two different fill densities were considered-1.2 g/cm 3 and 1.7 g/cm 3 and no measurable change was observed. Thus while the density of the fill is an unknown quantity, given that the fill density is visibly less than that of the low density clay whose density was measured (as above), the value of 1.2 g/cm 3 was used in final dose rate calculation. Thus, while the exact density of the Bêt Khallaf tomb fill can not be reconstructed accurately, simulation of two different scenarios illustrates that within the likely range of fill densities which could be present at Bêt Khallaf, variation within this range will not affect the modelled dose rate.
6. The vessel's proximity to the tomb walls, that is, its location within the tomb. Was it in the centre, corner, or along a side of the tomb?
A limited study of this was carried out at low resolution and no measurable changes were observed. Given the very low attenuation rates of gamma rays in air with low moisture content, this is to be expected.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: OSL Optically stimulated luminescence IMA Infinite matrix assumption DRAC Dose Rate and Age Calculator

Appendix A. Infinite Minute Matrix Assumption Plus Geometric Assumption
We wish to determine the radiation absorbed per unit time in an infinitesimal volume δV in the situation illustrated in Figure A1. In the case where Material 2 is identical to Material 1, the infinite matrix assumption (IMA) applies and the radiation absorbed is equal to E γ1 δV where E γ1 is the volumetric emissivity of Material 1. However, we may also write it as where r is the radius and Γ(r, φ) is the dose rate in δV that originates at radius r and solid angle φ. Now let us consider the case where Material 2 is not identical to Material 1. Let us first define Γ a 1 (r, φ) to be the dose rate from material a where only one collision occurs (equivalent to path γ 1 in Figure A1).
where f 1 and f 2 are the fractions of the total volume of the white sphere occupied by Materials 1 and 2 respectively. The right hand side of Equation (A2) is precisely the geometric approximation discussed in Section 1.2.2. What about if intermediate collisions cannot be neglected? In fact the direct path γ 1 typically accounts for only half the dose (Ref. [5]: 289). In this case the geometric assumption still holds approximately for Figure A1 provided the attenuation of gamma radiation is similar for Materials 1 and 2. This is likely to be the case if, for example, Materials 1 and 2 have similar bulk compositions, differing chiefly in their radioisotope profiles. If, on the other hand, the materials have different attenuation properties, then the geometric assumption is not expected to apply.
Lastly, we note that in the case where the sample does not lie at the interface between two materials, one may use approximations like those described in Ref. [5]: 291.

Appendix B. Determining Chemical Compositions
Appendix B.1. The Bulk Composition of the Clays The composition data for the Nile silt clay was provided by five bulk samples initially tested by a commercial laboratory, which in addition to the elemental analysis, provided a breakdown of the main components of the clay in percentages. With only five samples analysed it was not possible to know the specific chemical breakdown of each individual vessel from the Bêt Khallaf assemblage. Instead, based on these five samples, an average breakdown was determined (See Table A1). Table A2 gives the chemical composition breakdown of marl clay and Nile silt clay used in the DOSIVOX simulations (for more detailed discussion see Ref. [11]).  Table A2. Chemical breakdown of Nile clay, take from Table A1, and marl clay, taken from [24], used in the DOSIVOX simulations. The radioisotope concentrations in the clay of the vessels were also needed to calculate the internal dose rate, and thus were already available for inclusion in the DOSIVOX model. Owing to the samples being taken from museum specimens, the minimum extraction technique was used for nearly all sample collection [15]. The first ∼ 1 ⁄2 of the sample was analysed using ICP-MS. Special care was taken with sample collection and measurement owing to the small size of the sample [11,15].
Tables A3 and A4 give the measured concentrations of radioisotopes in the clay and (in the case of the Tomb of Djer) contents of the vessels from Bêt Khallaf and Tomb of Djer respectively.  For the Bêt Khallaf material, it was noticed during the sampling process that several vessels had material adhering to their interior surface (vessels X5458, X5458, X5463, X5467, and X4120), or lodged tightly up under the rim of the vessel (vessels X5458 and X5459). It was possible therefore to collect this material for ICP-MS analysis and use the radionuclide concentrations as a proxy for the entire assemblage (see Table A5). The ICP-MS values for this project that were used to determineḊ ext are presented in Table A5, with the average value being used. It was decided that an average value would be the best value to use forḊ ext calculation, even though the values varied slightly between the individual vessels yielding external sediment. This was because the long attenuation depth of the external gamma radiation means thatḊ ext is affected by a larger body of material than the sediment in direct contact with the vessel. In this case, an average value will be more representative than an individual measurement as it will compensate for the variation in composition of the surrounding sediment.