Monika Korte
email: monika@gfz-potsdam.de

Maxwell C Brown
email: maxwell@hi.is

S R Gunnarson

Andreas Nilsson
email: andreas.nilsson@geol.lu.se

Sanja Panovska
email: panovska@gfz-potsdam.de

Ingo Wardinski
email: ingo.wardinski@univ-nantes.fr

Catherine G Constable
email: cconstable@ucsd.edu

Refining Holocene geochronologies using palaeomagnetic records

Keywords: Geochronology, Palaeomagnetic dating, Archaeomagnetic dating, Palaeosecular Variation

The aperiodic nature of geomagnetic field variations, both in intensity and direction, can aid in dating archaeological artefacts, volcanic rocks, and sediment records that carry a palaeomagnetic signal. The success of palaeomagnetic dating relies upon our knowledge of past field variations at specific locations. Regional archaeo-and palaeomagnetic reference curves and predictions from global geomagnetic field models provide our best description of field variations through the Holocene. State-of-the-art palaeomagnetic laboratory practices and accurate independent age controls are prerequisites for deriving reliable reference curves and models from archaeological, volcanic, and sedimentary palaeomagnetic data. In this review paper we give an overview of these prerequisites and the available reference curves and models, discuss techniques for palaeomagnetic dating, and outline its limitations. In particular, palaeomagnetic dating on its own cannot give unique results, but rather serves to refine or confirm ages obtained by other methods. Owing

Introduction

Earth's magnetic field, also known as the geomagnetic field, is primarily generated in Earth's molten convecting metallic outer core. As fired archaeological materials and igneous rocks cool, and sediments are deposited at Earth's surface, assemblages of magnetic iron oxide minerals within them have the capacity to record the past geomagnetic field (known as the palaeomagnetic field). In addition, chemical sedimentary rocks, such as flowstones or speleothems, which form through the precipitation of minerals, have the potential to record the palaeomagnetic field and are used in an increasing number of studies (e.g. [START_REF] Lascu | Speleothem magnetism[END_REF][START_REF] Ponte | Speleothems as magnetic archives: paleosecular variation and a relative paleointensity record from a portuguese speleothem[END_REF][START_REF] Zanella | A 10,000 yr record of high-resolution paleosecular variation from a flowstone of Rio Martino Cave, Northwestern Alps, Italy[END_REF]. Although the age of the first magnetization recorded on Earth (∼4.2 Ga) is currently debated [START_REF] Tarduno | Antiquity of the south atlantic anomaly and evidence for top-down control on the geodynamo[END_REF][START_REF] Weiss | Secondary magnetic inclusions in detrital zircons from the Jack Hills, Western Australia, and implications for the origin of the geodynamo[END_REF] the Earth is likely to have had a magnetic field as long ago as 3.5 Ga [START_REF] Tarduno | Geodynamo, solar wind, and magnetopause 3.4 to 3.45 billion years ago[END_REF].

The geomagnetic field is a vector field and three convenient quantities to describe it are declination, inclination, and intensity (Fig. 1a). Declination is the angle of the vector eastward of geographic north in the horizontal plane.

Inclination is the angle between the vector and Earth's surface (positive downward). The field intensity (strength) is the magnitude of the vector.

Over its long history the geomagnetic field has varied through time, making its behaviour useful for geochronological applications. Core field variations occur on timescales of months to millions of years [START_REF] Constable | A paleomagnetic power spectrum[END_REF] and the field varies non-uniformly over the globe. The magnitude of these changes can vary significantly and there is a range of rates of change. Variations of the core magnetic field over months and longer, as observed in modern measurements, are known as secular variation. Variations in the palaeomagnetic field prior to direct measurements are named Figure 1: (a) Schematic representation of a magnetic field vector in a geographic coordinate system and definitions of declination, inclination, and intensity. (b) Geocentric dipole configuration for the geomagnetic field. In this example, geographic north is offset from geomagnetic north (the dipole axis is tilted) by a similar amount to today. Field lines are not shown in the source region of the geodynamo, the outer core. palaeosecular variation (PSV). The term in general describes variations on centennial timescales and longer, because higher frequency variations cannot be resolved from archaeological artefacts, volcanic rocks, and sediments that record the palaeomagnetic field. More extreme palaeomagnetic field changes are field reversals and excursions. Reversals are defined by a complete swap of the magnetic poles, i.e., north and south poles change position. Geomagnetic excursions also involve extreme field changes with some directions suggesting a partially or fully reversed polarity; however, unlike reversals, the field recovers to its original polarity rather than switching poles for a prolonged time period. Reduced intensity commonly accompanies the direction variations of both reversals and excursions (e.g., [START_REF] Roberts | Geomagnetic excursions: knowns and unknowns[END_REF][START_REF] Brown | Microwave palaeointensity results from the Matuyama-Brunhes geomagnetic field reversal[END_REF][START_REF] Channell | Reconciling astrochronological and 40 Ar/ 39 Ar ages for the Matuyama-Brunhes boundary and late Matuyama chron[END_REF][START_REF] Laj | Geomagnetic excursions[END_REF]. Reversals and excursions have been widely used as chronologically useful stratigraphic markers (e.g. [START_REF] Lowrie | Magnetostratigraphy[END_REF][START_REF] Florindo | Eocene-Oligocene magnetobiochronology of ODP Sites 689 and 690, Maud Rise, Weddell Sea, Antarctica[END_REF][START_REF] Mazaud | Geomagnetic-assisted stratigraphy and sea surface temperature changes in core MD94-103 (Southern Indian Ocean): possible implications for North-South climatic relationships around H4[END_REF][START_REF] Collins | Palaeomagnetic and biostratigraphic dating of marine sediments from the Scotia Sea, Antarctica: First identification of the Laschamp excursion in the Southern Ocean[END_REF]. In this paper, however, we focus on the Holocene (past ∼12 kyr) when such extreme geomagnetic variations did not occur.

PSV has been used to provide age constraints on archaeomagnetic materials, volcanic rocks, and sediments on timescales of hundreds to thousands of years. By determining the declination, inclination, and palaeointensity (the field strength in the past) recorded by materials with unknown ages, these parameters can be compared with known (independently dated) records of field change, and an age can be assigned. This method of dating is known as archaeomagnetic dating or palaeomagnetic dating and was originally proposed by [START_REF] Thellier | Sur l'aimantation des terres cuites et ses applications geophysiques[END_REF]. However, as often noted, palaeomagnetic field evolution cannot give unique age information because values of field intensity and directions recur over time (e.g., [START_REF] Aitken | Dating by archaeomagnetic and thermoluminescent methods[END_REF][START_REF] Thellier | Early research on the intensity of the ancient geomagnetic field[END_REF][START_REF] Clark | Developments in archaeomagnetic dating in Britain[END_REF], even though variations are not periodic. It has been suggested, e.g., by Emile Thellier (see [START_REF] Aitken | Dating by archaeomagnetic and thermoluminescent methods[END_REF] that "the term 'magnetic dating' should be used with reserve". The term magnetic age refinement is more precise because the method can only be applied to materials that can be allocated to a certain archaeological or geological epoch by other methods. Note, that on even longer timescales than considered here the term palaeomagnetic dating is also used when virtual geomagnetic pole positions from ancient rocks, combined to form apparent polar wander paths, or the geomagnetic polarity timescale (e.g. [START_REF] Cande | Revised calibration of the geomagnetic timescale for the late Cretaceous and Cenozoic[END_REF][START_REF] Ogg | Geomagnetic Polarity Time Scale[END_REF] are used for assigning or refining chronologies.

A complication in palaeomagnetic age refinement is that field variations are non-uniform across the globe (as illustrated in Section 2). Although the geomagnetic field is dominated by a simple dipole geometry, more complex fields, known as non-dipole fields (see Section 2), contribute to regional differences in field structure, which cannot be neglected if PSV is used for chronology. Geomagnetic field evolution can only be used to refine ages if a reliable PSV reference curve exists for the region of interest. Reference curves can be produced in several ways: (a) regional curves based on archaeomagnetic and volcanic data; (b) regional stacks of sediment records; (c) global stacks of sediment data, which may or may not incorporate constraints from absolute palaeointensity data; or (d) combining all available data in regional or global inverse magnetic field models. Note, that in the context of this paper and the timescales we are dealing with the term palaeomagnetic reference curve always refers to curves of continuous time variations in field direction and/or intensity. There are two main aspects that must be considered when developing a palaeomagnetic reference curve: (1) the reliability of the palaeomagnetic data determined in the laboratory, and (2) the accuracy and uncertainties of ages assigned to the data by, e.g., archaeological or radiometric methods. These aspects are discussed in Sections 3 and 4.

Archaeomagnetic dating curves have been used for age refinement since the 1960s and are now a common part of archaeological studies. In some countries, e.g., the UK, archaeomagnetic dating curves have been used by governmental heritage agencies as part of the procedure of age refinement for artefacts (see, [START_REF] Batt | Advances in archaeomagnetic dating in Britain: New data, new approaches and a new calibration curve[END_REF]. Some recent applications have been to refine the ages of archaeological structures (e.g., [START_REF] Carrancho | Fullvector archaeomagnetic dating of a Medieval limekiln at Pinilla Del Valle Site (Madrid, Spain)[END_REF][START_REF] Hammond | An archaeomagnetic study of a Roman bath in Southern France[END_REF][START_REF] Principe | Archaeomagnetic dating of Copper Age furnaces at Croce di Papa village and relations on Vesuvius and Phlegraean fields volcanic activity[END_REF], to place age constraints on human activities [START_REF] Peters | Archaeomagnetic dating of pyrotechnological contexts: a case study for copper smelting sites in the Central Timna Valley, Israel[END_REF], habitation of peoples in certain areas, e.g., through investigation of burial contexts (Goguitchaichvili et al., 2017b), and cultural practices (Goguitchaichvili et al., 2017a).

A growing number of reference curves and global geomagnetic field models have been published over recent years, with particular focus on either the most recent few millennia for which many archaeomagnetic data are available or extended to the entire Holocene where sediment records are dominant.

These curves and models offer new possibilities for palaeomagnetic chronology refinement, but caveats exist. In this paper, we review the possibilities and limitations of archaeo-and palaeomagnetic age refinement on centennial to millennial timescales during the Holocene.

Following an introduction to global and regional aspects of the geomagnetic field (Section 2), we briefly summarize the work-flow used to obtain high quality archaeo-and palaeomagnetic data and outline characteristics of different data types (Section 3). PSV can only be used for geochronological applications if an independent chronology of (regional) magnetic field evolution has first been established; Section 4 is dedicated to this requirement.

In Section 5 we give an overview of palaeomagnetic reference curves, stacks, and models and give examples of how they have been used to refine chronologies for different applications. Through an analysis of spatial and temporal field variability (Section 6), we provide guidelines on the limits of correlating palaeomagnetic data to reference records. Our conclusions include requirements for future improvements of geomagnetic field models and reference curves to serve as age refinement tools.

The geomagnetic field

In its simplest form the global geomagnetic field can be imagined as having a structure similar to a bar magnet placed at Earth's centre (Fig. 1b).

This structure is referred to as a dipole. Over geological timescales a dipolar field is considered a reasonable assumption and forms the basis of the geocentric axial dipole (GAD) approximation, where the dipole axis aligns with Earth's rotation axis. This configuration results in spatially simple patterns of declination, inclination, and intensity at Earth's surface (Fig. 2a). Under the GAD assumption, declination is zero across Earth's surface. Inclination, I, and field intensity, F , are related to latitude λ or colatitude θ (with θ = 90 • -λ) at Earth's radius r by tanI = 2tanλ

(1) and

F = µ 0 P 4πr 3 √ 1 + 3cos 2 θ, ( 2
)
where P is the dipole moment of the field and µ 0 the vacuum permeability. Inclination is positive in the northern hemisphere and negative in the southern hemisphere, with zero inclination along the equator, increasing in absolute magnitude toward the poles. Intensity is minimum at the equator and increases toward the poles.

In reality the dipole is rarely aligned with Earth's rotation axis. For example, applying the dipole assumption to today's field results in a deviation of the geomagnetic north pole (i.e., the axis of the approximated dipole field) by approximately 10 • from the geographic north pole, with its calculated location over Ellesmere Island, northern Canada. Although the surface field of the approximated tilted dipole is still spatially simple, the tilt of the dipole axis adds a small degree of complexity (Fig. 2b). Declination is non-zero, site dependent, and becomes larger in the polar regions. Inclination and intensity retain their latitudinal dependent structures; however, the magnetic equator undulates and is offset from the geographic equator.

Reducing the geomagnetic field to a (tilted) dipole is conceptually straightforward; however, the global field structure at Earth's surface has more complexity and spatial variation than can be explained by a dipole alone (Fig. 2c), as illustrated in maps derived from the International Geomagnetic Reference Field, 12 th generation (IGRF-12) for 2015 [START_REF] Thébault | International geomagnetic reference field: the 12th generation[END_REF]. The tilted dipole accounts for ∼93% of the present day geomagnetic main field power at Earth's surface (excluding lithospheric and external sources). The remainder after subtraction of the best-fit dipole or tilted dipole is referred to as the non-axial-dipole or non-dipole field, respectively. This consists of spatially smaller scale structures, which can be described mathematically, e.g., by a sum of increasingly small-scale spherical harmonic functions, as quadrupoles, octupoles, etc. [START_REF] Merrill | The use of magnetic field excursions in stratigraphy[END_REF] give visual examples of these structures. The field structure at Earth's surface can be thought of as a combination of fields with different geometries and magnitudes. Non-dipole field structures result in a global geomagnetic field that can vary significantly with location (Fig. 2c). This is important to consider when assessing field behaviour at globally disparate locations; a topic that will be revisited throughout this article. A clear example of how the non-dipole field can influence global field structures is the intensity minimum called the South Atlantic Anomaly that is visible in today's field (Fig. 2c). The influence of the non-dipole field results in an offset between the locations where the measured geomagnetic field is vertical (inclination = ±90

• ), known as the dip poles or magnetic poles, and the geomagnetic poles calculated assuming a geocentric dipole (Fig. 3). Geomagnetic poles will always be antipodal to each other, whereas non-dipole components result in north and south magnetic poles that are frequently not antipodal.

A comparison of palaeomagnetic results from different locations is not straightforward using directional or intensity data because the dipole contribution causes clear (mainly latitudinal) inclination and intensity differences, e.g., the field intensity at the equator today is about 30 µT, but is closer to 65 µT around the poles (Fig. 2a andb). The concepts of virtual dipole moment (VDM), virtual axial dipole moment (VADM), and virtual geomagnetic pole (VGP) can be used to eliminate geographic variations due to a dipole field from intensity and directional data, respectively. A VADM represents the dipole moment of a pure axial dipole field with the observed site intensity. It can be obtained from an intensity measurement F considering the site colatitude θ,

V ADM = 4πr 3 µ 0 √ 1 + 3cos 2 θ (3)
(see also eq. 2). A VDM represents the dipole moment of a tilted dipole, where information about the tilt comes from knowledge of the site inclination, mathematically expressed as site colatitude:

cotθ = tanI 2 (4)
in eq. 3 (see eq. 1). A VGP is the geomagnetic pole location if the observed field direction at a site was due to a pure tilted dipole field and can be calculated from the palaeomagnetic declination and inclination together with site latitude and longitude (see, e.g., [START_REF] Butler | Paleomagnetism: Magnetic domains to geologic terrains[END_REF]Merrill et al., 1996, for all equations). Conveniently, the programs di vgp.py and b vdm.py from the PmagPy collection of palaeomagnetic programs by [START_REF] Tauxe | PmagPy: Software package for paleomagnetic data analysis and a bridge to the Magnetics Information Consortium (MagIC) Database[END_REF] (Appendix B) can be used to transform directions to VGPs and intensity to VDM or VADM. If declination is not available, as is the case with azimuthally un-oriented sediment cores and partly oriented archaeological artefacts (e.g., bricks or potsherds where only the vertical, but not the horizontal orientation during firing can be estimated), the so-called inclination anomaly can be determined by subtracting the axial dipole inclination for the site latitude from inclination data. As illustrated in this section the GAD approximation, and transformations of the data to VDMs, VADMs, and VGPs can be extremely useful in removing the largest expected geographic field variations. However, we must note a significant limitation. Despite widespread application of the GAD approximation within the palaeomagnetic community, it is unclear over how long the field must be averaged or even whether this assumption is always valid. 10,000 years has often been considered long enough (e.g., [START_REF] Merrill | The Magnetic Field of the Earth: Palaeomagnetism, the Core, and the Deep Mantle[END_REF]. However, existing averages of both global Holocene data and field models depart systematically from the GAD average [START_REF] Constable | Persistent high paleosecular variation activity in Southern Hemisphere for at least 10 000 years[END_REF]. The most pronounced (and widely agreed upon) contribution is that of an average axial quadrupole, which contributes to stronger average fields in the northern hemisphere than in the south, but there are also indications of longitudinal field variability, with weaker average fields in the western Pacific than in the east. Additionally, there are signs of greater overall PSV activity in the southern hemisphere, and lower variability in the Pacific hemisphere [START_REF] Constable | Persistent high paleosecular variation activity in Southern Hemisphere for at least 10 000 years[END_REF]. Studies based on older lava flow data provide strong indications that such departures from GAD spatial structure extend to million year timescales [START_REF] Cromwell | PSV10: A global data set for 0-10 Ma time-averaged field and paleosecular variation studies[END_REF]. It is therefore important to consider that there remains no clear consensus on whether and on what timescale the GAD approximation can be considered adequate.

Archaeological material and volcanic rocks

Archaeological materials and volcanic rocks acquire a TRM on initial cooling or after reheating [START_REF] Dunlop | Rock magnetism: Fundamentals and frontiers[END_REF]. This can result in a strong and geologically stable remanence. This remanence is acquired instantly on geological timescales and allows recovery of well-defined snapshots of direction and intensity. The thermal nature of the remanence allows 'absolute' palaeointensity estimates to be made by replacing the original TRM with a laboratory TRM acquired by heating in the laboratory in a known field (see [START_REF] Tauxe | 13 -paleointensities[END_REF]. However, the archaeomagnetic and volcanic record is fragmented in time and each palaeomagnetic estimate requires an independent date. Lavas provide only spot readings in time and may lack any stratigraphic control, and there can often be hundreds to thousands of years between eruptions of lavas at the same or nearby locations.

Similarly, palaeomagnetic data from archaeological materials often come from archaeological sites that were active for only certain periods of time. Data density also depends on age. Data from archaeological materials are abundant over the past 2000 years, but decline rapidly for older ages, with few data prior to 10 ka (see Section 5, Fig. 4 and Brown et al., 2015b Volcanic and archaeological samples are recovered in various ways depending on material type. Providing that a sample is drilled or cut from an in-situ structure (e.g., a lava flow or kiln), the physical orientation of the sample may also be recorded, e.g., using a sun compass or a magnetic compass (e.g., Turner et al., 2015b;English Heritage, 2006). In that case the full vector palaeomagnetic signal might be determined, otherwise the sample can only be used for palaeointensity determination.

Sediments

In contrast to palaeomagnetic data from archaeological materials and volcanic rocks, the process by which sediments acquire a magnetization is not fully understood (see, e.g., [START_REF] Roberts | Magnetic paleointensity stratigraphy and high-resolution Quaternary geochronology: successes and future challenges[END_REF][START_REF] Tauxe | 13 -paleointensities[END_REF], despite considerable work to understand the issues that influence sediment remanence acquisition (e.g. [START_REF] Egli | Natural remanent magnetization acquisition in bioturbated sediment: General theory and implications for relative paleointensity reconstructions[END_REF][START_REF] Zhao | Microbially assisted recording of the Earth's magnetic field in sediment[END_REF][START_REF] Valet | Detrital magnetization of laboratoryredeposited sediments[END_REF][START_REF] Chen | Remanence acquisition efficiency in biogenic and detrital magnetite and recording of geomagnetic paleointensity[END_REF]. The overriding principle is that some fraction of magnetic particles aligns with the geomagnetic field after their last disruption event after deposition as sediment. The final alignment may occur near the sediment-water interface and the remanence is 'locked in' at this time (a detrital remanent magnetization) or it may continue during burial with the final lock-in occurring at some depth within the sediment (a postdepositional remanence).

Sediments are usually collected in gravity or piston cores, from which discrete cubes (typically 6 to 8 cm 3 ) or a long strip of sediment (known as a u-channel [START_REF] Tauxe | u channels: a new technique for paleomagnetic analysis of hydraulic piston cores[END_REF][START_REF] Nagy | New advances for paleomagnetic studies of sediment cores using u-channels[END_REF]) can be extracted by carefully pushing plastic samples into the split-half of a sediment core. In some locations it is possible to directly sample exposed sections of sediments.

On long timescales, sediment data are primarily from marine environments, whereas on shorter timescales such as during the Holocene sedimentary data sets are dominated by cores from lacustrine environments. The latter often have higher sedimentation rates that can resolve decadal to centennial scale magnetic field variations. The advantage that sediment palaeomagnetic records have over archaeomagnetic and volcanic data is that they are quasi-continuous and can span thousands to hundreds of thousands of years depending on the setting. An additional advantage is that the stratigraphic order of palaeomagnetic observations is known.

However, all sediment palaeomagnetic records are a smoothed representation of the field to some degree. This is a result of the remanence acquisition process [START_REF] Roberts | Why are geomagnetic excursions not always recorded in sediments? Constraints from post-depositional remanent magnetization lock-in modelling[END_REF][START_REF] Roberts | Magnetic paleointensity stratigraphy and high-resolution Quaternary geochronology: successes and future challenges[END_REF] and the method of measurement. Both of these issues depend on sedimentation rate: the lower the sedimentation rate, the greater the smoothing. In addition to natural processes that influence remanence acquisition, the method of sampling and measurement can also cause smoothing of the palaeomagnetic record. For a discrete cube sample the magnetometer measures the magnetization acquired over a period of time dependent on the sedimentation rate and the length of the specimen. For example, for a core with sedimentation rate of 20 cm/ka sampled with continuous discrete specimens at a 2 cm spacing, the averaging will be 100 years. For u-channel measurements the core is passed through the magnetometer and measurements are made at a set spacing [START_REF] Weeks | Improvements in long core measurement techniques: Applications in paleomagnetism and paleoceanography[END_REF][START_REF] Nagy | New advances for paleomagnetic studies of sediment cores using u-channels[END_REF]. As the sediment in the uchannel is a continuous sample and the magnetometer sensor has a response function of a few centimetre width, the sensor integrates the magnetization on either side of the point directly below the sensor. A highly smoothed record can result. Moreover, remanence acquisition in sediments may be delayed compared to sediment age as a result of post-depositional processes (see, e.g., [START_REF] Tauxe | 13 -paleointensities[END_REF][START_REF] Suganuma | 10 Be evidence for delayed acquisition of remanent magnetization in marine sediments: Implication for a new age for the Matuyama-Brunhes boundary[END_REF][START_REF] Roberts | Magnetic paleointensity stratigraphy and high-resolution Quaternary geochronology: successes and future challenges[END_REF][START_REF] Mellström | Post-depositional remanent magnetization lock-in depth in precisely dated varved sediments assessed by archaeomagnetic field models[END_REF][START_REF] Nilsson | Short-term magnetic field variations from the post-depositional remanence of lake sediments[END_REF]. This must be considered when determining age models for palaeomagnetic records.

Laboratory experiments and methods

The weak natural remanent magnetization (NRM) of palaeomagnetic samples is measured with highly sensitive cryogenic or spinner magnetometers (e.g., Turner et al., 2015b). Depending on the sample type it is common practice to use stepwise thermal or alternating field (AF) demagnetization to remove magnetic overprints associated with coring, transportation/storage, or secondary heating events. Demagnetization data are analysed using orthogonal projections (e.g. [START_REF] Zijderveld | demagnetisation of rocks: analysis of results[END_REF], where horizontal and vertical magnetization components are simultaneously projected and NRM components can be separated. To determine the direction of the "characteristic" remanent magnetization (ChRM), principal component analysis (PCA) is commonly used to find the best least-squares line fit through the demagnetization data [START_REF] Kirschvink | The least-squares line and plane and the analysis of palaeomagnetic data[END_REF]. However, in many early studies demagne-tization data were taken at a single AF or heating step once any viscous magnetization had been removed and the directional data appeared stable.

Such approaches are no longer used and detailed stepwise demagnetization is carried out to ensure that a reliable ChRM is obtained. The experimental error of a PCA line fit is represented by the maximum angular deviation (MAD), which is calculated from the variance of the data around the principal axis. Improvements to common practices for obtaining realistic uncertainty estimates from PCA have recently been proposed by Heslop and Roberts (2016a,b).

To assess the reliability of the palaeomagnetic signal and to determine the grain size and composition of magnetic components within a specimen (remanence carrying or otherwise), various rock magnetic analyses are commonly employed. In general, similar rock magnetic methods can be used for archaeological materials, volcanic rocks, and sediments. Thermomagnetic analysis (measurement of an induced magnetization or susceptibility across a range of temperatures), allows the Curie temperature (or temperatures)

of assemblages of magnetic grains to be determined and their mineral composition inferred. Through comparison of magnetization during heating and cooling, the propensity for thermal alteration can be assessed. This is important for understanding whether absolute palaeointensity experiments are likely to be successful [START_REF] Tauxe | 13 -paleointensities[END_REF]. Magnetic grain size can be determined by a number of methods, including hysteresis measurements (Jackson and Solheid, 2010;[START_REF] Paterson | Measuring, processing, and analyzing hysteresis data[END_REF], FORC diagrams [START_REF] Roberts | First-order reversal curve diagrams: A new tool for characterizing the magnetic properties of natural samples[END_REF][START_REF] Roberts | Understanding fine magnetic particle systems through use of first-order reversal curve diagrams[END_REF], AF demagnetization of NRM and anhysteretic remanent magnetization (ARM), and isothermal remanent magnetization (IRM)

acquisition. Over the past 20 years increasingly sophisticated methods and analysis techniques have been developed to investigate magnetic grain assemblages in rocks and sediments (e.g., [START_REF] Roberts | First-order reversal curve diagrams: A new tool for characterizing the magnetic properties of natural samples[END_REF][START_REF] Kruiver | Quantification of magnetic coercivity components by the analysis of acquisition curves of isothermal remanent magnetisation[END_REF][START_REF] Egli | Characterization of individual rock magnetic components by analysis of remanence cureves, 1. unmixing natural sediments[END_REF][START_REF] Lascu | Quantifying the concentration of ferrimagnetic particles in sediments using rock magnetic methods[END_REF][START_REF] Heslop | A method for unmixing magnetic hysteresis loops[END_REF][START_REF] Roberts | Understanding fine magnetic particle systems through use of first-order reversal curve diagrams[END_REF][START_REF] Heslop | Numerical strategies for magnetic mineral unmixing[END_REF]. The results of grain size analyses are frequently non-unique and this research area is currently at the forefront of rock magnetism [START_REF] Heslop | Numerical strategies for magnetic mineral unmixing[END_REF][START_REF] Roberts | A critical appraisal of the "Day" diagram[END_REF]. In addition, anisotropy of magnetic susceptibility (AMS) is useful for (1) assessing any disturbance in the sedimentary fabric and (2) as an initial test to assess whether there is a remanence anisotropy in archaeological materials. For archaeological materials, this can be followed up by anisotropy of ARM (AARM) and/or anisotropy of TRM (ATRM) measurements (e.g., [START_REF] Chauvin | Paleointensity of the geomagnetic field recovered on archaeomagnetic sites from France[END_REF]. This allows the remanence directions and intensity to be corrected for remanence anisotropy.

The methods used to determine palaeointensity from archaeological/volcanic materials and sediments are fundamentally different. The magnetic properties of sediments can be influenced by variations in magnetic grain concentration, grain size, mineralogy, lithology, and other environmental influences.

The bulk magnetic properties of sediments are commonly used to calculate relative palaeointensity (RPI) by normalising the NRM (e.g. NRM/ARM, NRM/SIRM) in an attempt to eliminate these influences and to extract the geomagnetic field signal [START_REF] King | A new rock magnetic approach to selecting sediments for geomagnetic paleointensity studies: application to paleointensity for the last 4000 years[END_REF][START_REF] Tauxe | Sedimentary records of relative paleointensity of the geomagnetic field: theory and practice[END_REF][START_REF] Roberts | Magnetic paleointensity stratigraphy and high-resolution Quaternary geochronology: successes and future challenges[END_REF][START_REF] Tauxe | 13 -paleointensities[END_REF]. However, no absolute field intensity information can be retrieved. Moreover, the normalization is not always successful, and attention to changes in the magnetic properties of the sediment is required.

Not all sediment records are suitable for RPI analysis [START_REF] King | A new rock magnetic approach to selecting sediments for geomagnetic paleointensity studies: application to paleointensity for the last 4000 years[END_REF][START_REF] Tauxe | Sedimentary records of relative paleointensity of the geomagnetic field: theory and practice[END_REF][START_REF] Roberts | Searching for single domain magnetite in the "pseudo-single-domain" sedimentary haystack: Implications of biogenic magnetite preservation for sediment magnetism and relative paleointensity determinations[END_REF][START_REF] Tauxe | 13 -paleointensities[END_REF].

Absolute palaeointensity can be obtained from materials that carry a TRM. Absolute palaeointensity determination, pioneered by [START_REF] Thellier | Sur l'intensité du champ magnétique terrestre dans le passé historique et géologique[END_REF], is based on successively replacing the NRM (assumed to be a TRM) at increasing temperature steps with a new magnetisation acquired in a known field (a laboratory-TRM). If the relationship between the NRM lost and the laboratory TRM gained at increasing temperature steps is linear, the gradient of the straight line fit to these data can be multiplied by the value of the applied laboratory field to give the ancient intensity. Variations on the original Thellier method (e.g., [START_REF] Coe | Paleo-intensities of the Earth's magnetic field determined from Tertiary and Quaternary rocks[END_REF][START_REF] Aitken | Determination of the intensity of the Earth's magnetic field during archaeological times: Reliability of the Thellier technique[END_REF][START_REF] Tauxe | A simplified statistical model for the geomagnetic field and the detection of shallow bias in paleomagnetic inclinations: was the ancient magnetic field dipolar? In: Timescales of the paleomagnetic field[END_REF] and alternative methods based upon different principles have been developed to overcome the influence of non-ideal remanence carriers (e.g., multi-domain grains) and sample alteration, which can result in experimental failure. Alternative methods include the Shaw method and derivatives [START_REF] Shaw | A new method of determining the magnitude of the palaeomagnetic field. Application to five historic lavas and five archaeological samples[END_REF][START_REF] Rolph | A new method of paleofield magnitude correction for thermally altered samples and its application to Lower Carboniferous lavas[END_REF][START_REF] Yamamoto | Palaeointensity study of the Hawaiian 1960 lava: implications for possible causes of erroneously high intensities[END_REF], the multispecimen method [START_REF] Dekkers | Reliable absolute palaeointensities independent of magnetic domain state[END_REF][START_REF] Fabian | Multiple-specimen absolute paleointensity determination: An optimal protocol including pTRM normalization, domain-state correction, and alteration test[END_REF], the microwave method (e.g., [START_REF] Walton | Application of ferrimagnetic resonance heating to palaeointensity determinations[END_REF][START_REF] Hill | Palaeointensity results for historic lavas from Mt Etna using microwave demagnetization/remagnetization in a modified Thellier-type experiment[END_REF], the Triaxe vibrating sample magnetometer method (Le [START_REF] Goff | A new three-axis vibrating sample magnetometer for continuous high-temperature magnetization measurements: applications to paleo-and archeo-intensity determinations[END_REF]) and the Wilson method [START_REF] Wilson | Paleomagnetism in Northen Ireland, Part I, The thermal demagnetization of natural magnetic moments in rocks[END_REF][START_REF] Muxworthy | Revisiting a domain-state independent method of palaeointensity determination[END_REF]. [START_REF] Perrin | Paleointensity determinations, magnetic domain structure, and selection criteria[END_REF][START_REF] Valet | Time variations in geomagnetic intensity[END_REF], [START_REF] Dunlop | Physical basis of the Thellier-Thellier and related paleointensity methods[END_REF], and [START_REF] Tauxe | 13 -paleointensities[END_REF] provide overviews of experimental procedures and the issues surrounding absolute palaeointensity determinations. In addition, improvements have been made in our understanding of selecting the most favourable materials and grain sizes for obtaining successful palaeointensity results (e.g., [START_REF] Carvallo | Increasing the efficiency of paleointensity analyses by selection of samples using first-order reversal curve diagrams[END_REF][START_REF] Valet | A selective procedure for absolute paleointensity in lava flows[END_REF][START_REF] Cromwell | Paleointensity estimates from historic and modern Hawaiian lava flows using glassy basalt as a primary source material[END_REF][START_REF] Paterson | Bulk magnetic domain stability controls paleointensity fidelity[END_REF] and in methods for analysing and determining the reliability of experimental data [START_REF] Paterson | On improving the selection of Thellier-type paleointensity data[END_REF].

Additional data treatments

To be able to use palaeomagnetic data, in particular sediment records, to refine age estimates or to constrain geomagnetic field models it is often necessary to introduce additional data treatments of raw directions and RPI values. In the following sections we briefly describe the most common types of data treatment. It is important to note that all of these treatments rely on assumptions, whether it be in sampling, laboratory procedures or geomagnetic field behaviour. Therefore, the raw data should always be reported.

Core orientations and relative declination calibration

Sediment cores are usually not oriented in the horizontal plane (with respect to geographic north) and declination data are, therefore, typically reported as relative values (e.g., Turner et al., 2015b). In cases where several cores are retrieved in a sequence it is common to rotate the lower section into alignment with the upper section so that the declination data produce a continuous sequence across the core break. Such inter-core adjustments are often subjective, or even not possible [START_REF] Geiss | A Holocene-Late Pleistocene geomagnetic inclination record from Grandfather Lake, SW Alaska[END_REF], but can be greatly improved if cores are retrieved with a slight depth overlap to the previous section.

Core rotation while it penetrates the sediments is likely a common but frequently ignored problem [START_REF] Snowball | Geomagnetic field intensity changes in Sweden between 9000 and 450 cal BP: extending the record of "archaeomagnetic jerks" by means of lake sediments and the pseudo-Thellier technique[END_REF], which can depend on the type of corer used (Turner et al., 2015b). Rotations, which affect declination data, can be detected by comparing overlapping cores [START_REF] Stanton | Assessing the reliability of Holocene relative palaeointensity estimates: a case study from Swedish varved lake sediments[END_REF] and can potentially be corrected by detrending the data, e.g. assuming a constant rotation rate [START_REF] Ali | Holocene palaeomagnetic secular variation at Lake Biwa, central Japan[END_REF].

For direct comparison or inclusion in field models, it is necessary to reorient records to absolute declination. Often records are oriented so that they have zero mean over their whole time interval and are referred to as relative declinations [START_REF] Turner | Lake sediment record of the geomagnetic secular variation in Britain during Holocene times[END_REF]. Technically this is the same as assuming an average GAD field. However, the validity of the GAD assumption in general is unclear [START_REF] Cromwell | PSV10: A global data set for 0-10 Ma time-averaged field and paleosecular variation studies[END_REF] and it cannot easily be justified for comparatively short time intervals of a couple of thousand years or less (see section 2). Orientation by absolute declinations from oriented materials or models are preferable. However, care should be taken when such orientations are based on comparisons between historical field observations and palaeomagnetic data from the uppermost, typically slushy, sediments that often carry a less stable palaeomagnetic signal and are more prone to physical deformation during coring (e.g., [START_REF] Nourgaliev | Very high resolution paleosecular variation record for the last ∼ 1200 years from the Aral Sea[END_REF].

Sediment cores that do not penetrate the sediment vertically, which affects both inclination and declination, are another frequently encountered problem [START_REF] Constable | Holocene geomagnetic secular variation records from north-eastern Australian lake sediments[END_REF][START_REF] Stoner | A paleomagnetic approach toward refining Holocene radiocarbon-based chronologies: paleoceanographic records from the north Iceland (MD99-2269) and east Greenland (MD99-2322) margins[END_REF][START_REF] Stanton | Assessing the reliability of Holocene relative palaeointensity estimates: a case study from Swedish varved lake sediments[END_REF].

If sediments contain horizontal interbeds or laminations it may be possible to measure the angle directly on the core [START_REF] Stanton | Assessing the reliability of Holocene relative palaeointensity estimates: a case study from Swedish varved lake sediments[END_REF]. However, in more or less homogeneous sediments the issue is only identified due to inconsistencies between palaeomagnetic data from parallel cores [START_REF] Turner | A 5000 year geomagnetic palaeosecular variation record from western Canada[END_REF]Turner et al., 2015a). To correct both inclination and declination data for tilted core penetration it is necessary to rotate the data with an appropriate pole and rotation angle, e.g., on a unit sphere. In cases where the tilt is unknown, [START_REF] Denham | Numerical correlation of recent paleomagnetic records in two Lake Tahoe cores[END_REF] devised a "palaeomagnetic pattern matching"

technique to find the pole and angle that maximizes pair-wise correlation of two unit vectors. This technique has been used to recover consistent palaeomagnetic signals in multiple cores that initially had up to 20-30 • deviations in mean inclination [START_REF] Constable | Holocene geomagnetic secular variation records from north-eastern Australian lake sediments[END_REF][START_REF] Turner | A 5000 year geomagnetic palaeosecular variation record from western Canada[END_REF].

RPI scaling

Calibration of sediment RPI to absolute values has been attempted in several ways. RPI series reaching into modern times can be compared with direct magnetic field observations or models thereof, but the time overlap is generally too short for robust calibration. More commonly, direct comparisons of RPI to absolute values determined on archaeomagnetic samples or volcanic rocks from nearby locations have been made (e.g., [START_REF] Constable | Holocene geomagnetic secular variation records from north-eastern Australian lake sediments[END_REF][START_REF] Constable | Palaeointensity in the pelagic realm: marine sediment data compared with archaeomagnetic and lake sediment records[END_REF][START_REF] Donadini | Geomagnetic field for 0-3 ka: 1. New data sets for global modeling[END_REF]. Over longer timescales than the Holocene (several hundreds of kyr), RPI stacks have been calibrated by comparison to globally averaged volcanic VADM results (Guyodo [START_REF] Guyodo | Global changes in intensity of the Earth's magnetic field during the past 800 kyr[END_REF][START_REF] Valet | Geomagnetic dipole strength and reversal rate over the past two million years[END_REF] and individual RPI records have been calibrated by making assumptions about the non-dipole field strength during a reversal [START_REF] Constable | Towards absolute calibration of sedimentary paleointensity records[END_REF]. Neither of these two methods is applicable to calibrating RPI records on Holocene timescales.

An alternative approach is to calibrate Holocene RPI records by comparison with regional or global geomagnetic field models, where scaling factors are estimated as the median of the ratio between the whole time series (Korte and Constable, 2006). In this case, it must be carefully considered how well constrained the model is for the desired region. For determining global models, available RPI records are mostly scaled prior to modelling, based on intensities from nearby regions (e.g., [START_REF] Donadini | Geomagnetic field for 0-3 ka: 1. New data sets for global modeling[END_REF] or on a previous model that is constrained by absolute intensities (Korte and Constable, 2011;[START_REF] Licht | Ensembles of low degree archeomagnetic field models for the past three millennia[END_REF][START_REF] Nilsson | Reconstructing Holocene geomagnetic field variation: new methods, models and implications[END_REF]. In some cases the calibration is refined iteratively in the modelling process (Korte and Constable, 2011).

An attempt to co-estimate RPI calibration factors during the inversion with inclusion of absolute intensities without prior scaling [START_REF] Panovska | Limitations in paleomagnetic data and modelling techniques and their impact on Holocene geomagnetic field models[END_REF] demonstrated that so far the number of available absolute absolute intensities spanning the Holocene is insufficient for robust calibration. [START_REF] Nilsson | Reconstructing Holocene geomagnetic field variation: new methods, models and implications[END_REF] noted that some Holocene sediment records contain jumps due to depositional environment changes through time when calibrated by model curves. Different scaling factors are appropriate for different parts of the sequence in such cases.

U-channel data instrument response

Narrow access pass-through magnetometers enable rapid and dense measurements, typically every cm, of u-channel samples. However, the shape of the response function of the pick-up coils, rather than the measurement density, dictates the actual resolution of each measurement, with typical halfpower widths in the range of 5-8 cm [START_REF] Weeks | Improvements in long core measurement techniques: Applications in paleomagnetism and paleoceanography[END_REF][START_REF] Nagy | New advances for paleomagnetic studies of sediment cores using u-channels[END_REF]Jackson et al., 2010;[START_REF] Oda | Toward robust deconvolution of pass-through paleomagnetic measurements: new tool to estimate magnetometer sensor response and laser interferometry of sample positioning accuracy[END_REF]. In addition, response functions for the transverse and axial measurement axes are typically not the same, which leads to both a smoothed and distorted measurement of the NRM [START_REF] Roberts | High-resolution magnetic analysis of sediment cores: strengths, limitations and strategies for maximizing the value of long-core magnetic data[END_REF]. Several techniques have been developed to correct for these issues by deconvolving the data according to known system responses [START_REF] Oda | Deconvolution of long-core paleomagnetic data of ocean drilling program by Akaike's Bayesian Information Criterion minimization[END_REF]Jackson et al., 2010;[START_REF] Oda | Deconvolution of continuous paleomagnetic data from pass-through magnetometer: a new algorithm to restore geomagnetic and environmental information based on realistic optimization[END_REF], with the recently developed software UDECON facilitating implementation of such algorithms [START_REF] Xuan | UDECON: deconvolution optimization software for restoring high-resolution records from pass-through paleomagnetic measurements[END_REF]. However, it remains difficult to attain a continuous record because of measurement end effects so that data tend to be truncated in the upper and lowermost 5 cm of u-channel samples.

Assigning uncertainties

All palaeomagnetic and most chronological data have uncertainties. Uncertainty estimates provide important information because they account for the accuracy and precision of measurements and ensure a proper weighting of data in reference curves and models. Many factors can influence the accuracy of palaeomagnetic and archaeomagnetic data, and there is no straightforward way to assign uncertainties. Often no error estimates are reported for palaeomagnetic data, particularly for sediments where there often is only one measurement per depth.

Archaeomagnetic and volcanic directions are typically reported as a mean obtained using Fisher statistics [START_REF] Fisher | Dispersion on a sphere[END_REF] from multiple independently oriented samples from the same time horizon (e.g., from multiple samples from a hearth or a lava flow). Uncertainties are then based on dispersion of the directions and are reported as a radius of confidence around the mean direction at the significance 95% level [START_REF] Fisher | Dispersion on a sphere[END_REF], also known as α 95 (see [START_REF] Butler | Paleomagnetism: Magnetic domains to geologic terrains[END_REF], for a more detailed explanation). α 95 can be converted to standard deviation errors of declination and inclination (see, e.g., [START_REF] Donadini | Geomagnetic field for 0-3 ka: 1. New data sets for global modeling[END_REF]. Similarly, uncertainties on mean absolute intensities are mostly given by the standard deviation from multiple samples.

For a single palaeomagnetic sediment core, two recent studies proposed techniques to transform the experimental error from principal component analysis of directional data into meaningful errors on individual results. The approach of [START_REF] Khokhlov | Principal component analysis of palaeomagnetic directions: converting a maximum angular deviation (MAD) into an α95 angle[END_REF] µT intensity uncertainties were used to construct the Holocene CALSx models (e.g., Korte and Constable, 2011;Korte et al., 2011). [START_REF] Suttie | Consistent treatment of errors in archaeointensity implies rapid decay of the dipole prior to 1840[END_REF] found that systematic errors appear to be an important factor in archaeointensity data. This, in turn, means that using percentages of intensities as uncertainties produces biased estimates.

New uncertainty estimates for sediment palaeomagnetic records have been obtained by investigating scatter around fits using a robust smoothing spline technique [START_REF] Panovska | Spline analysis of Holocene sediment magnetic records: Uncertainty estimates for field modeling[END_REF]

Independent chronology

A palaeo-or archaeomagnetic time series with a good chronology that has been obtained independently from magnetic field variations is a prerequisite for a reference curve intended to be used for palaeomagnetic age refinement.

Several methods provide independent ages for Holocene archaeological ma- [START_REF] Telford | The intercept method is a poor estimate of a calibrated radiocarbon age[END_REF][START_REF] Michczyński | Is it possible to find a good point estimate of a calibrated radiocarbon date?[END_REF][START_REF] Blaauw | Testing the timing of radiocarbon-dated events between proxy archives[END_REF][START_REF] Scott | Radiocarbon dating: sources of error[END_REF]. It is, therefore, common to report an age range at the 2 σ significance level in archaeological and volcanic studies, although more sophisticated modelling approaches are used in some archaeological studies.

Radiocarbon dating is the most frequently used method for dating Holocene sediments (see review by [START_REF] Zimmerman | Lacustrine Environments (14C)[END_REF]. When creating a sediment record the sequence is initially assigned a depth scale, which is subsequently transformed into time. It is not possible to measure the age at every depth in the sediment core, so ages are generally determined for some age-depth tie-points across the record. To transfer depth to time to produce a time series for a core, some form of interpolation is required between calibrated radiocarbon dates. Many approaches haven been taken to construct age-depth models (see [START_REF] Björck | 14 C chronostratigraphic techniques in paleolimnology[END_REF][START_REF] Telford | The intercept method is a poor estimate of a calibrated radiocarbon age[END_REF][START_REF] Parnell | A flexible approach to assessing synchroneity of past events using Bayesian reconstructions of sedimentation history[END_REF][START_REF] Bronk Ramsey | Deposition models for chronological records[END_REF][START_REF] Blaauw | Methods and code for 'classical' age-modelling of radiocarbon sequences[END_REF][START_REF] Blaauw | Flexible paleoclimate age-depth models using an autoregressive gamma process[END_REF]. The simplest approach, often found in older publications, takes point estimates of age (mean or median calibrated or uncalibrated age) and applies piecewise linear interpolation. Such models imply that sediments are deposited at the same rate between ages and that sedimentation rate changes are abrupt rather than gradual. This approach is also insensitive to hiatuses.

As point estimates are used, this method assumes that the age at a specific depth is precisely known and has no uncertainty. These models are not statistically meaningful. As noted by [START_REF] Björck | 14 C chronostratigraphic techniques in paleolimnology[END_REF], if one were to simply consider ±1σ uncertainty on a set of ages, then only 68.3% of the ages would lie within their ±σ uncertainties. Uncertainties are non-zero, so there is no statistically reasonable case for interpolating between points. Applying linear or polynomial regression fitting generally assumes that uncertainties are normally distributed; this is not the case for calibrated radiocarbon ages with a probabilistic distribution with many options for interpolation between them. Building a depth-age model using Bayesian statistics is the best strategy for resolving these issues (e.g. Bronk Ramsey, 1995; [START_REF] Blaauw | Flexible paleoclimate age-depth models using an autoregressive gamma process[END_REF]. This approach makes multiple depth-age models that fit a set of initial parameters and calculates uncertainty bounds across all depths.

Varying approaches for treating calibrated ages and the methods used to construct age-depth models may result in differences among palaeomagnetic time series purely from the methods employed. This must be considered when TL dating can be applied to fired archaeological materials such as pottery (see [START_REF] Aitken | Science-based dating in archaeology[END_REF][START_REF] Duller | Luminescence Dating[END_REF]. However, TL dating coupled with archaeomagnetic results has only been attempted in a limited number of studies (23 studies report TL dates in GEOMAGIA50). The precision of TL dating is less than for radiocarbon dating, with uncertainties typically in the range of 5-10% at ±1σ (68%) confidence limits [START_REF] Aitken | Science-based dating in archaeology[END_REF][START_REF] Bailiff | Luminescence, pottery and bricks[END_REF]. and volcanic data and also related age information. It is strongly encouraged to upload newly published magnetic data including dating information.

Information

PSV as a method for chronology refinement

In principle, any independently dated PSV record can be used to refine the age or chronology of a material with palaeomagnetic information. However, a palaeomagnetically refined chronology can be no more precise than the independent age information that went into the PSV reference curve used. Increasing numbers of independently dated archaeo-and palaeomagnetic results have led to development of reference curves, stacks, and global and regional models that can increase the reliability of PSV features and their ages through cross-validation of consistent signals from different sources.

Magnetic field variations may seem to recur in individual components, although PSV is not periodic. The resulting non-uniqueness in palaeomagnetic dating is reduced if two field components or ideally the full vector field information can be used for comparison with a reference. The time interval and region of interest also play a role in available options and the obtainable degree of palaeomagnetic age or chronology refinement.

Reference curves and models

Development of archaeomagnetic regional composite curves began in the mid-20 th century to increase the precision of archaeological chronologies.

Such work mainly focussed on particular countries, e.g., Japan (e.g., [START_REF] Watanabe | Secular variation in the direction of geomagnetism as the standard scale for geomagnetochronology in Japan[END_REF][START_REF] Hirooka | Archaeomagnetic study for the past 2,000 years in southwest Japan[END_REF], Britain (e.g., Aitken, 1958;[START_REF] Clark | Developments in archaeomagnetic dating in Britain[END_REF],

France [START_REF] Thellier | Sur la direction du champ magnétique terrestre, en France, durant les deux derniers millénaires[END_REF], Bulgaria [START_REF] Kovacheva | Summarized results of the archaeomagnetic investigation of the geomagnetic field variation for the last 8000 yr in south-eastern Europe[END_REF], and the United States of America (e.g., [START_REF] Dubois | Secular variation in southwestern United States as suggested by archaeomagnetic studies[END_REF][START_REF] Eighmy | Archaeomagnetic dating in the American Southwest[END_REF][START_REF] Sternberg | Archaeomagnetic secular variation of direction and paleointensity in the American Southwest[END_REF]. Independently dated archaeomagnetic results from spatially distributed locations (mostly spanning several hundred to a few thousand km, see Section 6 for a discussion of spatial correlation of PSV features) must be reduced to a central reference location to eliminate systematic differences owing to the dominant dipole field geometry. This can be achieved by relocating directions via their VGP (see, e.g., [START_REF] Shuey | Geographic correction of archaeomagnetic data[END_REF][START_REF] Noel | A method for correcting geographically separated remanence directions for the purpose of archaeomagnetic dating[END_REF]) and intensity through the VDM or VADM (see, e.g., [START_REF] Daly | An updated and homogeneous world secular variation data base. 1. smoothing of the archeomagnetic results[END_REF]. These assumptions do not allow non-dipole field differences to be accounted for. Construction of a regional reference curve generally requires a compromise between the size of the considered area and temporal data coverage. It takes years to decades of careful field and laboratory measurements to obtain enough individual archaeomagnetic and/or volcanic results to construct good reference curves. Many reference curves have been updated over time and new ones have been developed as additional archaeomagnetic results became available.

Recently compiled archaeomagnetic reference curves exist for several European countries, parts of Asia, and western North America (see Constable and Korte, 2015, for an overview), which roughly reflects the currently available data distribution (cf. Fig. 4b). They range from two (e.g., [START_REF] Márton | Two thousand years of geomagnetic field direction over central Europe revealed by indirect measurements[END_REF][START_REF] Yoshihara | Secular variation of the geomagnetic field intensity during the past 2000 years in Japan[END_REF] to eight [START_REF] Kovacheva | Geomagnetic field variations as determined from Bulgarian archeomagnetic data. Part II: the last 8000 years[END_REF] millennia in length and in general cover either declination and inclination (e.g., [START_REF] Gallet | Three millennia of directional variation of the Earth's magnetic field in western Europe as revealed by archaeological artefacts[END_REF][START_REF] Schnepp | Archaeomagnetic secular variation in Germany during the past 2500 years[END_REF][START_REF] Zananiri | Archaeomagnetic secular variation in the UK during the past 4000 years and its application to archaeomagnetic dating[END_REF][START_REF] Hagstrum | Achaeomagnetic dating in western North America: an updated reference curve based on paleomagnetic and archeomagnetic data sets[END_REF], or field intensity (e.g., [START_REF] Marco | New archaeointensity results from archaeological sites and variation of the geomagnetic field intensity for the last 7 millennia in greece[END_REF][START_REF] Yoshihara | Secular variation of the geomagnetic field intensity during the past 2000 years in Japan[END_REF], and only rarely all three components (e.g., [START_REF] Kovacheva | Geomagnetic field variations as determined from Bulgarian archeomagnetic data. Part II: the last 8000 years[END_REF].

Archaeomagnetic reference curves vary in the method of their construction and presentation. Directional data may be shown as individual time series (e.g., [START_REF] Watanabe | Secular variation in the direction of geomagnetism as the standard scale for geomagnetochronology in Japan[END_REF][START_REF] Kovacheva | Summarized results of the archaeomagnetic investigation of the geomagnetic field variation for the last 8000 yr in south-eastern Europe[END_REF], VGPs (e.g., [START_REF] Sternberg | Techniques for constructing secular variation curves and for interpreting archaeomagnetic dates[END_REF][START_REF] Labelle | Additional archaeomagnetic data on the south-west USA master geomagnetic pole curve[END_REF][START_REF] Lengyel | The pre-AD 585 extension of the U.S. Southwest archaeomagnetic reference curve[END_REF], or as Bauer plots [START_REF] Bauer | On the secular motion of a free magnetic needle[END_REF], where inclination is plotted against declination (e.g., [START_REF] Thellier | Sur la direction du champ magnétique terrestre, en France, durant les deux derniers millénaires[END_REF][START_REF] Clark | Developments in archaeomagnetic dating in Britain[END_REF][START_REF] Gallet | Three millennia of directional variation of the Earth's magnetic field in western Europe as revealed by archaeological artefacts[END_REF][START_REF] Hagstrum | Achaeomagnetic dating in western North America: an updated reference curve based on paleomagnetic and archeomagnetic data sets[END_REF].

Traditional simple interpolation of individual components or directional data in the Bauer plot with various methods with or without consideration of data uncertainties has now mainly been replaced by more sophisticated methods. [START_REF] Gallet | Three millennia of directional variation of the Earth's magnetic field in western Europe as revealed by archaeological artefacts[END_REF] proposed a bivariate extension to standard circularly symmetric [START_REF] Fisher | Dispersion on a sphere[END_REF] statistics to construct directional reference curves and to consider the temporal data density by adjustable moving windows. They applied the method to obtain a reference curve for France. Thébault and Gallet (2010) used a bootstrap method to derive an ensemble of curves providing a probability distribution for a Middle Eastern reference curve. A Bayesian method introduced by Lanos ( 2004) and [START_REF] Lanos | Hierarchical modelling of archaeomagnetic data and curve estimation by moving average technique[END_REF] accounts for both magnetic and age uncertainties and the uneven temporal data distribution and produces consistent curves for one to three field components with 95% confidence limits. The method has been widely applied, particularly to many recent reference curves for European countries (e.g., [START_REF] Schnepp | Archaeomagnetic secular variation in Germany during the past 2500 years[END_REF][START_REF] Tema | Archaeomagnetism in Italy: a compilation of data including new results and a preliminary Italian reference curve[END_REF][START_REF] Zananiri | Archaeomagnetic secular variation in the UK during the past 4000 years and its application to archaeomagnetic dating[END_REF][START_REF] Márton | Two thousand years of geomagnetic field direction over central Europe revealed by indirect measurements[END_REF]. Hellio et al. ( 2014) applied a Bayesian method using a time correlation function derived from present day geomagnetic field time spectra as prior information for intensity data from Syria and directional data from France. Their probability density reference curves recover more rapid variations than previous curves, which tended to be smoother owing to data and age uncertainties.

Well-dated sediment records can in principle be used individually as references for palaeomagnetic chronology refinement (e.g., Ólafsdóttir et al., 2013). If several records exist for a region they can be stacked to improve the reliability of the resulting reference curve. Examples are the British master curve [START_REF] Thompson | British geomagnetic master curve 10,000-0 yr BP for dating European sediments[END_REF][START_REF] Turner | Detransformation of the British geomagnetic secular variation record for Holocene times[END_REF] and the Fennoscandian directional and intensity stacks [START_REF] Snowball | FEN-NOSTACK and FENNOPRIS: Varve dated Holocene palaeomagnetic secular variation and relative palaeointensity stacks for Fennoscandia[END_REF]. On longer timescales (beyond the Holocene), global or regional intensity stacks have been widely used for chronological purposes (e.g., [START_REF] Laj | High resolution global paleointensity stack since 75 kyr (glopis-75) calibrated to absolute values[END_REF][START_REF] Valet | Geomagnetic dipole strength and reversal rate over the past two million years[END_REF][START_REF] Channell | Stacking paleointensity and oxygen isotope data for the last 1.5 Myr (PISO-1500)[END_REF][START_REF] Ziegler | PADM2M: a penalized maximum likelihood model of the 0-2 Ma palaeomagnetic axial dipole moment[END_REF]. [START_REF] Roberts | Magnetic paleointensity stratigraphy and high-resolution Quaternary geochronology: successes and future challenges[END_REF] reviewed possibilities and limitations of stacked records in detail.

Several global VADM reconstructions exist that average all available archaeomagnetic and volcanic intensity data at the time the curves were constructed. They span several millennia as illustrated in Fig. 5 [START_REF] Yang | Variations in the geomagnetic dipole moment over the last 12000 years[END_REF][START_REF] Genevey | ArcheoInt: An upgraded compilation of geomagnetic field intensity data for the past ten millennia and its application to the recovery of the past dipole moment[END_REF][START_REF] Knudsen | Variations in the geomagnetic dipole moment during the Holocene and the past 50 kyr[END_REF][START_REF] Valet | Secular variation of the geomagnetic dipole during the past 2000 years[END_REF]. Even when they are split into regions [START_REF] Genevey | ArcheoInt: An upgraded compilation of geomagnetic field intensity data for the past ten millennia and its application to the recovery of the past dipole moment[END_REF] or complemented by dipole tilt information [START_REF] Valet | Secular variation of the geomagnetic dipole during the past 2000 years[END_REF] they are less suitable for chronological purposes than reference curves or more detailed global models on these timescales, as they neglect regional non-dipole field differences. Moreover, their temporal resolution tends to be lower due to time averaging to minimize influences from non-dipole field contributions. The same is also true for past geomagnetic pole evolution reconstructions obtained by averaging VGPs from sediment records [START_REF] Nilsson | Holocene geocentric dipole tilt model constrained by sedimentary paleomagnetic data[END_REF][START_REF] Nilsson | Millennial scale cyclicity in the geodynamo inferred from a dipole tilt reconstruction[END_REF].

When using reference curves for dating one must consider the distance from the curve location (section 6) and it might be necessary to relocate directional data through the VGP equation to the coordinates of the dating curve and to calculate VADMs prior to attempting age refinement. These steps are not necessary when using curves predicted from models. Temporally continuous global magnetic field models that include the dipole and nondipole field in principle provide secular variation curves for any location on Earth, making them a convenient tool for magnetic chronology refinement. A growing number of models spanning recent centuries to millennia (Fig. 5) and built from historical, archaeo-, and palaeomagnetic data have been developed

over the last two decades. The models can be global or regional, in the latter case following similar methods to those described below for global models, either using different kinds of spatial basis functions [START_REF] Pavón-Carrasco | A regional archeomagnetic model for Europe for the last 3000 years, SCHA.DIF.3K: Applications to archeomagnetic dating[END_REF], 2010) or implemented as global models but limited to regional validity [START_REF] Lodge | Towards a new approach to archaeomagnetic dating in Europe using geomagnetic field modelling[END_REF][START_REF] Batt | Advances in archaeomagnetic dating in Britain: New data, new approaches and a new calibration curve[END_REF]). Both regional and global models can provide curves of all field components for any location within their region of validity and do not require assumptions about field geometry to compare directions and intensity.

The longest model that relies solely on direct magnetic field observations is guf m1 [START_REF] Jackson | Four centuries of geomagnetic secular variation from historical records[END_REF], which spans 400 years from 1590 to 1990. It is based on historical and modern magnetic field measurements and is dominated by declination information over the first two centuries (Fig. 4d). Many historical observations originate from shipboard measurements for navigational purposes, which provide a relatively good global data coverage (Fig. 4a). Absolute field intensity measurements only exist from the 1830s onward. An intensity scaling factor is required to build a global reconstruction based purely upon directional information [START_REF] Hulot | Uniqueness of mainly dipolar magnetic fields recovered from directional data[END_REF].

This factor was prescribed in the model by the strength of the axial dipole field contribution, linearly extrapolating its observed decrease from 1840 to 1990 to earlier times.

Early archaeo-and palaeomagnetic models [START_REF] Ohno | Spherical harmonic analysis of paleomagnetic secular variation curves[END_REF][START_REF] Hongre | An analysis of the geomagnetic field over the past 2000 years[END_REF][START_REF] Constable | Global geomagnetic field models for the past 3000 years: transient or permanent flux lobes? Phil[END_REF] provided only snapshots in time.

The next generation of models were temporally continuous, starting with a Continuous model from Archaeomagnetic and Lake Sediments for 3 kyr (CALS3K.1) [START_REF] Korte | Continuous global geomagnetic field models for the past 3000 years[END_REF]. Important improvements to subsequent models came from expanded archaeomagnetic, volcanic, and sedimentary data sets. Two types of models exist in terms of underlying data.

All three data types are used in the CALSxk (e.g., [START_REF] Korte | Continuous global geomagnetic field models for the past 3000 years[END_REF]Korte et al., 2011;[START_REF] Constable | Persistent high paleosecular variation activity in Southern Hemisphere for at least 10 000 years[END_REF], pfm9k [START_REF] Nilsson | Reconstructing Holocene geomagnetic field variation: new methods, models and implications[END_REF], and HFM [START_REF] Panovska | Limitations in paleomagnetic data and modelling techniques and their impact on Holocene geomagnetic field models[END_REF][START_REF] Constable | Persistent high paleosecular variation activity in Southern Hemisphere for at least 10 000 years[END_REF] families of models, the ASD FM and ASDI FM [START_REF] Licht | Ensembles of low degree archeomagnetic field models for the past three millennia[END_REF] models and the COV-LAKE model [START_REF] Hellio | Time-correlation-based regression of the geomagnetic field from archeological and sediment records[END_REF]. Only archaeomagnetic and volcanic data are used in the ARCHxk [START_REF] Donadini | Geomagnetic field for 0-3 ka: 1. New data sets for global modeling[END_REF][START_REF] Constable | Persistent high paleosecular variation activity in Southern Hemisphere for at least 10 000 years[END_REF] and S(C)HA (Pavón-Carrasco et al., 2008b,a, 2014b) families and the A FM [START_REF] Licht | Ensembles of low degree archeomagnetic field models for the past three millennia[END_REF], AmR [START_REF] Sanchez | Modelling the archaeomagnetic field under spatial constraints from dynamo simulations: a resolution analysis[END_REF] and COV-ARCH [START_REF] Hellio | Time-correlation-based regression of the geomagnetic field from archeological and sediment records[END_REF] models. The reasons for ignoring sedimentary magnetic field information in the latter models are to avoid issues of signal smoothing, lock-in delay, calibration of relative intensity, and correlation of temporal uncertainties.

However, this comes at the cost of a significantly reduced data set (Fig. 4), which leads to models that are poorly constrained by data in the southern hemisphere and weakly constrained globally for times prior to 1000 BCE.

Methodologically, nearly all existing continuous global models are based on expansions in spherical harmonic functions in space and on cubic splines in time (see, e.g. [START_REF] Jackson | Four centuries of geomagnetic secular variation from historical records[END_REF][START_REF] Korte | Continuous global geomagnetic field models for the past 3000 years[END_REF]. They are obtained by an inversion using, e.g., a regularized least-squares fit. Although models can be produced from archaeo-and palaeomagnetic data that fit the data within their uncertainty estimates, these tend to be overly complicated and dominated by smaller scale features than seen even in the present day field. The opposite is expected given the more limited field information compared with present day observations from geomagnetic observatories and magnetic satellite missions. A regularization is, therefore, implemented that limits the model resolution. The simplest form is a truncation of the spherical harmonic expansion at low degrees. However, because smaller scale signals present in the data cannot be accounted for in such a model, they might partly be mapped into the low degree model coefficients and distort them.

Most models now follow the philosophy of accommodating more structure than expected to be resolved in the basis function expansion, with an additional regularization constraint that minimizes some field quantity and leads to a smooth, simple model (see, e.g., [START_REF] Jackson | Four centuries of geomagnetic secular variation from historical records[END_REF]. A regularization factor trades off fit to the data against smoothness, to produce in this case a more dipole-dominated model [START_REF] Donadini | Geomagnetic field for 0-3 ka: 1. New data sets for global modeling[END_REF]. The factor is chosen such that the models do not contain smaller-scale structure than the present-day field [START_REF] Lodge | Towards a new approach to archaeomagnetic dating in Europe using geomagnetic field modelling[END_REF][START_REF] Donadini | Geomagnetic field for 0-3 ka: 1. New data sets for global modeling[END_REF], and an additional temporal regularization limits time variability. Differences mainly exist in data selection and data weighting by their uncertainties, outlier rejection, calibration of relative palaeointensities, orientation of relative declination, measure of misfit of the model to the data, maximum spherical harmonic degree, and method and strength of regularization.

Uneven and sometimes sparse data coverage, data and age uncertainties, and the consequently required regularization, lead to limited temporal resolution of all (presently available) archaeo-/palaeomagnetic models and reference curves. None of them fit large amplitude rapid field variations such as archaeomagnetic jerks [START_REF] Gallet | On the possible occurrance of 'archaeomagnetic jerks' in the geomagnetic field over the past three millennia[END_REF] or intensity spikes [START_REF] Ben-Yosef | Geomagnetic intensity spike recorded in high resolution slag deposit in Southern Jordan[END_REF][START_REF] Shaar | Geomagnetic field intensity: How high can it get? How fast can it change? Constraints from Iron Age copper slag[END_REF]. An example are the very well documented Levantine intensity spikes around 980 BCE and 740 BCE. The data suggest field variations that are more rapid than found in the present-day and historical field, and their geophysical origin is not understood [START_REF] Livermore | Core-flow constraints on extreme archeomagnetic intensity changes[END_REF][START_REF] Davies | Geomagnetic spikes on the core-mantle boundary[END_REF][START_REF] Korte | Archeomagnetic intensity spikes: global or regional geomagnetic field features? Front[END_REF]. The latter study shows that such localized rapid changes in archaeomagnetic data sets are generally not fit by standard interpolation methods used for models and reference curves.

In general, the largest differences between models come from data selection and weighting, and not from modelling strategy [START_REF] Panovska | Limitations in paleomagnetic data and modelling techniques and their impact on Holocene geomagnetic field models[END_REF][START_REF] Sanchez | Modelling the archaeomagnetic field under spatial constraints from dynamo simulations: a resolution analysis[END_REF], as can be seen in Fig. 6 where predictions from several models are compared for regions with different data coverage. Europe has dense archaeomagnetic and good sediment data coverage (Fig. 4). Africa, on the other hand, was nearly devoid of data when the latest published models were derived. Note that the situation has started to improve since then (see [START_REF] Korte | Global geomagnetic field reconstructions from centuries to excursions[END_REF], for a recent overview over newly published southern hemisphere data). Reasonable coverage with both archaeomagnetic and sediment data exists for Asia, whereas for southern South America a few sediment records but no archaeomagnetic/volcanic data (except for a handful of intensity data further north) were available when the models were built. Good agreement exists among all models for the past 3 kyr in Europe and with one exception the agreement is reasonably good in Asia for the same time interval. For other regions and times there mostly is general agreement in multi-centennial to millennial trends, but notable differences exist in variation details, particularly in declination and inclination. Not surprisingly, these differences are more pronounced when fewer data exist to constrain regional non-dipole field evolution. Models that do not include sediment data in several cases deviate strongly from general trends in models that incorporate all data types and are often less consistent with each other (e.g., declination in Europe around 1500 BCE), owing to the reduced amount of data constraining them.

Providing uncertainty estimates on model outputs is not straightforward. Model uncertainties are commonly estimated from ensembles of models obtained by re-sampling the data within their uncertainties (e.g., [START_REF] Donadini | Geomagnetic field for 0-3 ka: 1. New data sets for global modeling[END_REF][START_REF] Licht | Ensembles of low degree archeomagnetic field models for the past three millennia[END_REF][START_REF] Nilsson | Reconstructing Holocene geomagnetic field variation: new methods, models and implications[END_REF]Pavón-Carrasco et al., 2014b).

These model uncertainties, however, underestimate the uncertainty originating from a lack of information. They are mostly notably smaller than differences between individual models where and when models appear dis-

Examples of palaeomagnetic age and chronology refinement

Palaeomagnetic dating is applied either to refine individual ages mainly of archaeological or volcanic material, or to refine the chronology of a time series of magnetic results, mainly from sediment cores. These applications pose different challenges that may limit the accuracy of the dating improvement. Ambiguities in age due to recurring magnetic field values play an

Archaeomagnetic or volcanic age refinement

The easiest way to perform archaeomagnetic dating is by visual comparison of the measured magnetic directions or intensity to a reference curve.

This method was widely used in earlier studies, but is not satisfactory given the uncertainties on both palaeomagnetic measurements and calibration curves [START_REF] Batt | The British archaeomagnetic calibration curve: an objective treatment[END_REF]. Statistical methods that can take uncertainties into account and provide information on the accuracy and reliability of the determined age have been suggested, e.g., by [START_REF] Lanos | Bayesian inference of calibration curves: application to archaeomagnetism[END_REF] and Le [START_REF] Goff | On archeomagnetic secular variation curves and archeomagnetic dating[END_REF], (see also [START_REF] Mcintosh | An introduction ot archeomagnetic dating[END_REF]. The method suggested by Le Goff et al. ( 2002) is applicable to directional data and uses a bivariate extension of [START_REF] Fisher | Dispersion on a sphere[END_REF] statistics together with a test for whether two Fisherian distributions share a common mean direction [START_REF] Mcfadden | Classification of the reversal test in palaeomagnetism[END_REF].

The method of [START_REF] Lanos | Bayesian inference of calibration curves: application to archaeomagnetism[END_REF] provides probability density functions for possible dates for each palaeomagnetic component, which can be combined to give the most probable age with estimated uncertainty. [START_REF] Lanos | Bayesian inference of calibration curves: application to archaeomagnetism[END_REF] AD [START_REF] Kovacheva | Archaeomagnetic dating of archaeological sites from Switzerland and Bulgaria[END_REF]. Two additional intervals about 1000 years older and exceeding the 95% confidence limit illustrate that PSV can only refine chronological information in combination with other information. This example also demonstrates how knowledge of more than one field component can narrow the most probable ages (compare the combined probability density function with the individual ones in the middle row of Fig. 7). In this example, the most probable ages are mainly constrained by directional data while intensity alone gives inconclusive results.

Archaeomagnetic dating has been applied to better understand the archaeological history of sites all over the world with a broad range of ages.

Several recent examples demonstrate how the method can provide valuable age constraints, but on its own is often ambiguous or inconclusive. For instance, [START_REF] Tema | Combined archaeomagnetic and thermoluminescence study of a brick kiln excavated at Fontanetto Po (Vercelli, Northern Italy)[END_REF] and [START_REF] Casas | Archaeometric dating of two limekilns in an industrial heritage site in Calders (Catalonia, NE Spain)[END_REF] determined the last use of kilns in historical times (18 th to 20 th century) from inclination and intensity results from a site in Italy, and directional results from a site in Spain, respectively, using archaeomagnetic and historical models and confirmed by TL dating. However, the two magnetic components available to Tema et al.

(2013) gave two potential time intervals and hence larger uncertainty than the TL dating. A more conclusive application of the method was using full vector magnetic results from two kilns in Greece [START_REF] Tema | Dating of ancient kilns: A combined archaeomagnetic and thermoluminescence analysis applied to a brick workshop at Kato Achaia, Greece[END_REF]. They provided a rather narrow range of probable ages between 100 BCE and 100 CE using a model for the past 3 kyr, which is comparable to and in good agreement with dates determined by TL dating. Further examples from European sites where directional data helped to confirm or refine age ranges based on archaeological chronologies or radiocarbon include, e.g., kilns from Roman, Medieval and Modern times from Spain and Belgium [START_REF] Ech-Chakrouni | Constraints of archaeomagnetic dating and field intensity determinations in three ancient tile kilns in Belgium[END_REF][START_REF] Casas | Using archaeomagnetism to improve the dating of three sites in Catalonia (NE Spain)[END_REF] or burnt cave sediments from the Bronze age in Spain [START_REF] Carrancho | Fullvector archaeomagnetic dating of a Medieval limekiln at Pinilla Del Valle Site (Madrid, Spain)[END_REF].

As discussed above, past magnetic field variations are less well established for other parts of the world, so it is not surprising that attempts to apply archaeomagnetic dating can be less successful there. [START_REF] Lengyel | Archaeomagnetic research in the Andean highlands[END_REF] attempted to use a directional comparison with the historical gufm1 model to date furnaces, hearths, and cooking structures in the Bolivian Andes. In (2015) is another example for an area where palaeomagnetic age refinement remains difficult. Based on results from historically dated Chilean lava flows, they concluded that presently available models seem unsuitable for dating in South America due to the lack of regional data, and that archaeomagnetic dating there should be restricted to the last 3 centuries where magnetic field variation is known from direct historical observations.

Sedimentary chronology refinement

Sedimentary chronology refinement relies on correlating directional and/or relative palaeointensity variations to a reference curve. In many studies this is achieved purely by visual comparison of distinctive directional or intensity features. For example, it has been common to alphabetically annotate (with either Roman or Greek letters) key declination and inclination swings observed in reference records (e.g., [START_REF] Creer | Geomagnetic variations for the interval 7000-25[END_REF][START_REF] Creer | Geomagnetic secular variation in Greece through the last 6000 years obtained from lake sediment studies[END_REF][START_REF] Thompson | British geomagnetic master curve 10,000-0 yr BP for dating European sediments[END_REF][START_REF] Turner | Detransformation of the British geomagnetic secular variation record for Holocene times[END_REF][START_REF] Ojala | Testing the fidelity of sediment chronology: comparison of varve and paleomagnetic results from Holocene lake sediments from central Finland[END_REF][START_REF] Snowball | FEN-NOSTACK and FENNOPRIS: Varve dated Holocene palaeomagnetic secular variation and relative palaeointensity stacks for Fennoscandia[END_REF]. Similar features in the record to be dated have then been linked to these annotated features. More recently, computer algorithms have been developed to automatically match similar features within sediment records, e.g., the Match algorithm by [START_REF] Lisiecki | Application of dynamic programming to the correlation of paleoclimate records[END_REF]. An alternative approach is to incorporate PSV data directly into Bayesian age-depth models to provide additional temporal constraints. This method presented by [START_REF] Nilsson | Short-term magnetic field variations from the post-depositional remanence of lake sediments[END_REF] has the advantage that it uses all available data rather than a few identified PSV features and also takes data uncertainties into account. Beside providing a seamless way to combine PSV age refinement with conventional dating techniques (e.g., radiocarbon) the method also offers the possibility to correct for potential 'lock-in' delay in the palaeomagnetic signal.

The success of palaeomagnetic age refinement using sedimentary records varies, as indicated by the following recent examples. High-resolution palaeomagnetic records with robust geochronological controls from near the record in question can be used reliably (e.g., [START_REF] Ólafsdóttir | Synchronizing Holocene lacustrine and marine sediment records using paleomagnetic secular variation[END_REF][START_REF] Roza | Paleomagnetic secular variation and environmental magnetism of Holocene-age sediments from Tulare Lake, CA[END_REF]. Such correlation can be used to determine sedimentation rate variations to interpret environmental changes. [START_REF] Stoner | A paleomagnetic approach toward refining Holocene radiocarbon-based chronologies: paleoceanographic records from the north Iceland (MD99-2269) and east Greenland (MD99-2322) margins[END_REF] demonstrated how this approach can be used to transfer available radiocarbon dates between two records to refine age-depth models for both sites. Barletta et al. son with the Scandinavian directional stack FENNOSTACK [START_REF] Snowball | FEN-NOSTACK and FENNOPRIS: Varve dated Holocene palaeomagnetic secular variation and relative palaeointensity stacks for Fennoscandia[END_REF]. They used them to infer sedimentation rates discussed in the context of crustal uplift. They concluded that radiocarbon ages also obtained for the cores were systematically too old by up to 2500 years. [START_REF] Lougheed | Using an independent geochronology based on palaeomagnetic secular variation (PSV) and atmospheric Pb deposition to date Baltic Sea sediments and infer 14 C reservoir age[END_REF] similarly used PSV data compared with FENNOSTACK to infer down-core variations in the radiocarbon reservoir age of dated foraminifera in two sediment Baltic Sea cores.

Spatial and Temporal Correlation

The spatial extent over which correlations of PSV features for age refinements are reasonable is limited by the complexity of regional field morphology due to the non-dipole field and will, therefore, vary for different locations and time periods (e.g., [START_REF] Shuey | Geographic correction of archaeomagnetic data[END_REF][START_REF] Noel | A method for correcting geographically separated remanence directions for the purpose of archaeomagnetic dating[END_REF][START_REF] Casas | Distribution analysis of errors due to relocation of geomagnetic data using the Conversion via Pole (CVP) method: implications on archaeomagnetic data[END_REF]. In general, correlation breaks down with increasing distance. In areas (or times) with steep geomagnetic field gradients, e.g., the present-day South Atlantic, large spatial variations of different field components could lead to large chronological errors if care is not taken. Models of the surface field (Sections 5.1 and 6.1) and surface geomagnetic fields generated by numerical dynamo simulations (Section 6.2) can improve our understanding of spatial and temporal field variations. This, in turn, will aid determination of when PSV correlation is valid for chronological refinement.

Correlation based on Holocene geomagnetic field models

We can use Holocene geomagnetic field models to determine the distance over which two PSV records can be expected to reflect the same geomagnetic field behaviour on centennial to millennial timescales (e.g., when differences are smaller than data uncertainties). To get a first idea of regional field dif- differences to the central record increase with distance in all components.

The effect is stronger for directional than for intensity data and is strongest for declinations at high latitude locations close to the magnetic pole.

We have extended this analysis by considering global variation of palaeomagnetic time series away from a central location (again Windermere) for the 10 kyr duration of CALS10k.2. To do this, we produced palaeomagnetic time series of declination, inclination, and intensity over the globe with a grid spacing of 2 • in longitude and 1.5 • in latitude away from the central location. For each field component we calculated correlation coefficients between the Windermere time series and the individual time series of the global grid.

Correlation coefficients range from -1 to 1, where 1 is a perfect correlation, and intensity we used a Pearson correlation coefficient and for declination we used a circular correlation coefficient [START_REF] Fisher | Correlation coefficients for random variables on a unit sphere or hypersphere[END_REF]. To show global variations of the correlation coefficients we plot them as global contour maps, where correlations are shown at the coordinates of the gridded time series that were compared with the Windermere record (Fig. 11).

The non-symmetric correlation patterns seen in Fig. 11 so careful review of model data density for the selected location is advised (e.g., Fig. 1 in Korte et al., 2011). In some cases, model predictions may still be preferred over reference records or stacks. At high latitudes, for example, where traditional dating methods may be problematic due to environmental conditions (e.g., [START_REF] Wolfe | Geochronology of high latitude lake sediments[END_REF][START_REF] Lisé-Pronovost | Paleomagnetic constraints on the holocene stratigraphy of the arctic alaskan margin[END_REF][START_REF] St-Onge | Paleomagnetism near the north magnetic pole: A unique vantage point for understanding the dynamics of the geomagnetic field and its secular variations[END_REF], many studies have opted for direct comparisons of PSV records with CALSxk model predictions [START_REF] Lisé-Pronovost | Paleomagnetic constraints on the holocene stratigraphy of the arctic alaskan margin[END_REF][START_REF] Barletta | Dating of Holocene western Canadian Arctic sediments by matching paleomagnetic secular variation to a geomagnetic field model[END_REF][START_REF] Ledu | Holocene paleoceanography of the northwest passage, Canadian Arctic Archipelago[END_REF] tions, thus, may represent one end-member in correlations: smoothed records that are (presumably) largely dipole-dominated and that may be correlated over large distances at the cost of lost local and non-dipole field variations.

Insights from dynamo simulations

To overcome limitations with Holocene field models and to investigate correlation length based on a higher resolution long-term model we use an Earth-like dynamo simulation (Case A), to provide synthetic time series.

We compute similar spatial correlation maps as in Section 6.1. Details of this numerical simulation are given in Appendix A. In Figure 12 we show temporal evolution of the dipole tilt generated by Case A. The total run time of this numerical simulation is equivalent to about 2.4 Ma, i.e., 60 magnetic diffusion times (MDT). The simulation indicates numerous field excursions and at least 5 field reversals. The reversal rate matches approximately the reversal rate of Earth's magnetic field over the past 150 million years as given by the [START_REF] Cande | Revised calibration of the geomagnetic timescale for the late Cretaceous and Cenozoic[END_REF] polarity timescale. We consider temporal and spatial variations of this simulation to be broadly similar and probably slightly more variable than those for Earth. For our correlation analyses, we choose a period of stable field polarity that lasts ∼263 kyr. Results from this analysis can be thought of as the other end member in correlation when compared with the CALS10k.2 correlation analysis: conservative estimates on spatial coherence and resulting maximum expected errors for dating due to the dynamo simulation having higher spatial complexity than the presentday geomagnetic field. An additional advantage of using a dynamo simulation instead of a Holocene field model is that we can generate a large number of plausible geomagnetic field time series to get a statistical measure of spatial correlations of field components at different locations on Earth. There are no external forcing mechanisms in the Case A dynamo, so we do not expect persistent non-zonal (longitudinal) field variations nor persistent differences between northern and southern hemisphere field variations. It is, therefore, sufficient to investigate spatial correlation differences at different latitudes.

The normal-SV section of the dynamo simulation run was divided into 43 time-series of 100 simulated times steps each, roughly equivalent to 6000 years. For each time-series we followed the same steps to produce spatial correlation maps as in Fig. 11, but instead of using the location of Windermere the maps were calculated for latitudes of 0, 15, 30, 45, 60, 75, and 90 • N and longitude 0 • E and the corresponding coordinates in the opposite hemisphere (i.e., 0, 15, 30, 45, 60, 75, and 90 • S, and 180 • E). The distance between locations at the same latitude but in different hemispheres is large enough that the spatial correlation maps can be considered more or less independent from each other while the lack of persistent hemispherical asymmetries in the dynamo simulation means that they are statistically indistinguishable.

In addition to declination, inclination, and intensity we calculated correlation maps based on VGP latitude. Combining correlation maps from both hemispheres, we obtain N = 86 samples for each latitude and field component. To visualize the results (Fig. 13), we define a parameter R95 as the 5 th percentile of the 86 correlation coefficients R in each map grid cell (i.e., the expected lower limit of R with 95% confidence). Intensity time series at mid latitudes have the largest distances over which high correlations are obtained, which confirms the strong dipole influence in this component and its preference over directional information if no nearby reference curve is available. Differences in latitudinal correlation lengths could be explained by differences in the amount of vector field information contained in the intensity component of a dipole dominated field: at mid latitudes the scalar intensity is influenced strongly by all three orthogonal field components (north, east, and vertical), at the equator the influence is intermediate (mostly north and east, with a small vertical component), and it is lowest at the poles (mostly from the vertical component).

Some studies focused on high latitude PSV records seem to suggest that the high amplitude of inclination and declination features recorded close to the geomagnetic poles make PSV correlation especially robust in these locations (e.g., Ólafsdóttir et al., 2013). However, this is only true over short distances, as applied in that paper. From the figures presented here, both inclination and declination time-series quickly become spatially very incoherent at high latitudes. This is due to the proximity to the north and south magnetic poles, which introduces site-specific variations in the individual magnetic field direction. In extreme cases, e.g. with the pole moving back and forward between two nearby sites, this could even lead to anti-correlated signals (in this case in inclination). VGP transformation of inclination and declination (see Section 2) is based on a dipole field approximation and, therefore, removes a large part of the latitudinal field dependence. This probably explains why spatial correlation maps for VGP are more or less independent of latitude (Fig. 13). Using VGP latitudes also circumvents problems associated with inclination and declination at high latitudes because it measures the geomagnetic pole position instead of the site-specific magnetic field direction toward the pole. We elaborate on this explanation based on two examples in the following section.

VGP correlations and dipolar field variations

VGP latitude can be correlated over greater distances than declination and inclination alone (Fig. 13). VGP correlation distances are similar regardless of site location, whereas for declination and inclination correlation distances taper off significantly at high latitudes. The latter is due to shorter distances between the geomagnetic pole and observation site at high latitudes and we demonstrate here with two examples how large a role variable orientation of the tilted dipole plays even without non-dipole field contributions.

In Fig. 14, we show a hypothetical example of how a tilted, but otherwise purely dipolar, field is observed depending on the observation site. We simulated a change of the tilted dipole field in time by simply varying the geomagnetic pole (GP) longitude between 180 • and -180 • , and by keeping latitude constant. In this purely dipolar case the GP path is identical to a VGP path determined at any location on Earth. From this VGP path we calculated corresponding declination and inclination variations at two locations 45 • apart. There are two important outcomes of this simple exercise:

(1) moving the GP around a line of constant latitude varies both inclination and declination at an observation site;

(2) spatially separated sites have longitudinally offset inclination and declination variations, which correspond to temporal offsets. Inclination and declination variations between sites (for constant VGP latitude) occur because the dipole axis is either tilted toward or away from the observation site. This is apparent if we consider the inclination of site 1. When the VGP is directly north of site 1 (at 0 • longitude), the inclination is maximum as the dipole is tilted toward the site. Conversely, when the VGP is at 180 • longitude, the pole is at its furthest point away from site 1 and the inclination is a minimum. At these two points the declination is zero. When the VGP is directly north of site 1, it is north-east of site 2, and therefore gives a smaller inclination at site 2 than at site 1. As the apparent merely because of the relationship between the observation site and the geomagnetic pole location.

We can apply this concept to a real data set and assess the influence of temporal changes when changing the observation site. In Fig. 15a Results of this simple exercise demonstrate that using declination and inclination for long-distance correlation can be ambiguous (particularly at high latitudes) and great care must be taken with this approach. Alternatively, declination and inclination can be transformed to VGP latitude and longitude; all the declination and inclination variation in Fig. 15c andd then collapses to the same VGP latitude (Fig. 15b) and longitude curves.

Note that VGP variation used in this example will likely contain notable contributions from the non-dipole field and that the real tilted dipole field contribution probably generally varies less on comparable timescales, so that the effect from pure dipole tilt might be overestimated here. Regional non-dipole field differences that add to observed regional declination and inclination differences and tend to vary faster than the dipole are, on the other hand, neglected in this example. The observed directional differences in Fig. 15c andd could be of realistic magnitude, however, for the combined effect from dipole axis variations and non-dipole fields. If these curves came from different locations we could expect both effects to contribute to the difference. Non-dipole contributions would, therefore, mean that the conversion to VGPs would not yield identical curves as in this example. However, the difference between VGP records would be reduced in comparison to using declination and inclination records on their own. For VGP correlations to be successful, treatment of declination data from sediment cores needs to be carefully considered. To mitigate the lack of azimuthal orientation, declination data are commonly rotated so that the mean of the record is zero. Such a correction is required prior to making the VGP calculation; however, the duration of the records being compared must be sufficiently long so that the field averages to a zonal field, which might never be the case. If not, an incorrect zero-mean declination rotation will propagate into the VGP calculation and result in differences in VGP latitude and longitude calculated from both records. [START_REF] Lisiecki | Application of dynamic programming to the correlation of paleoclimate records[END_REF], URL given in Appendix B) with a high degree of flexibility, meaning that relatively short-term maxima and minima are matched, as is often done when selecting tie points for visual correlation.

The beginning and end points of the records are fixed in this case, which simulates the likely scenario that a palaeomagnetic time-series would have at least two well-established age control points between which to match PSV variations. It can be seen that the high correlation results from similarity of the general trends of the time-series. On short timescales, some features in the curves get offset in time by the matching. Consequently, matching individual components introduces, on average, centennial-scale errors when directional (declination and inclination) records are matched in this way. Intensity correlations, by comparison, generally have average offsets of >500 years. This difference is likely due to the higher-frequency PSV features observed in declination and inclination records which make for a good match on Holocene or shorter timescales. Intensity correlations might be more useful, however, on longer timescales and over greater distances due to the largely Assuming that examples from the smoothed inverse model CALS10k.2 (Fig. 16) might underestimate age errors, we also use predictions from the numerical simulation. We compare time series from numerical dynamo simulation runs (N = 86 for inclination, declination, and intensity; see Section 6.2) at 45 • latitude in both hemispheres for the equivalent of 6,000 years, with correlation coefficients of R≈0.8. In order to estimate the maximum expected errors for visual matching, only the beginning points of the time series are fixed, and endpoints of the target records were allowed to vary. Centennialscale mean offsets for all three components (inclination, declination, and intensity) are shown in Fig. 17. For both data-based models and dynamo simulations, maximum age offsets for all three components are millennial in scale (for inclination and declination) and greater than millennial-scale (for intensity). Although the dynamo simulation runs do not span the Holocene and do not specifically represent Holocene field behaviour, general agreement in error magnitude for field components gives a general precision of PSV tie points on Holocene timescales. However, note that components in this example are matched individually. In reality, all available component records from one sediment core must agree in the depth-age model after palaeomagnetic age refinement.

The precision of palaeomagnetic age refinement can be improved by increasing the number of independently and accurately dated tie points [START_REF] Mcmillan | Limitations in correlation of regional relative geomagnetic paleointensity[END_REF], or by matching two or three field components at once. Furthermore, the appropriateness of fit between a palaeomagnetic time series and a reference PSV curve depends on what is being investigated and the timescale in question, and no single method applies to all PSV records.

Estimates given here assume no error in the target and reference records. In reality, errors introduced only by correlation of PSV features which are assumed to originate from a similar field structure at a location should be propagated through from the original, independently dated reference chronology and/or field model, and may be higher or lower depending on the precision of the reference curve chronology and the resolution of the PSV curve being refined.

Conclusions and outlook

Palaeomagnetic directional and intensity variations can aid in dating archaeological artefacts, volcanic rocks and sediments that record the palaeomagnetic field. In this paper, we provided an overview of how palaeomagnetic age refinement has been applied in archaeological, volcanic, and sedimentary contexts from historical to Holocene times and discussed the potential of the method and its limitations. Prerequisites for age refinement are reliably determined palaeomagnetic data from the material to be dated and the existence of a high-resolution regional palaeomagnetic reference record that is robustly constrained in time.

Magnetic field variations are globally well known only over historical times, with declination and inclination well constrained for the past ≈300 to 400 years. Direct absolute intensity measurements, however, have only been made since 1840. Global models including all available historical data can be used to predict magnetic field variations reliably for any location on Earth for chronological purposes over these time intervals. Further back in time, our knowledge of the geomagnetic field relies on archaeo-or palaeomagnetic data. Large enough numbers of data with good independent age controls must be available for an area in order to establish a reliable regional reference curve. This is not the case globally. Spherical harmonic field models based on archaeo-and/or palaeomagnetic data give curves of direction and intensity for any location on Earth, but their reliability varies regionally depending on data coverage. Millennial and longer timescale reference curves and models do not resolve the full field variability. Variations in the model curves are smoothed and depend on measurement and dating uncertainties of the underlying data and additionally, for sediments, on sedimentation rate combined with sample size and lock-in time. All these factors will have an influence on the success of palaeomagnetic age refinement and must be considered when choosing an appropriate reference curve.

Use of a reference curve from a global model has the advantage that it can be obtained directly for the site of interest, whereas the nearest regional reference curve might have been developed for a location several tens to thousands of kilometres away. In this case, care must be taken to ensure that the material to be dated and the reference curve recorded the same geomagnetic signal, and that non-dipole field contributions can be neglected.

Strict guidelines cannot be given, but correlation analyses suggest that age errors can easily reach several centuries if only one field component is used and the distance between site location and reference curve exceeds a few hundred km. We recommend that data and reference curve are converted to VGP for directional data and VDM or VADM for intensity when their locations are not the same. Re-location to the site of the reference curve via the dipole field assumption, as is frequently done in archaeomagnetism, has the same effect. In general directions are more suitable for correlation of fast variations over short distances, whereas for larger distances (over ∼ several hundred km) and longer periods, intensities are preferable.

For several regions on Earth, particularly large parts of the Southern Hemisphere, our knowledge of Holocene geomagnetic field variations is not yet detailed enough for reliable chronological application. Progress will come from new data that add further information about past geomagnetic field variations. In particular, strong community efforts to produce new Southern Hemisphere and low latitude data will provide reference curves for additional regions and improve global models. Attempts to incorporate uncertainties into global models will give more realistic model uncertainty estimates that improve our understanding of the accuracy achievable in archaeo-or palaeomagnetic age refinement for individual regions. Important prerequisites are state-of-the-art archaeo-and palaeomagnetic laboratory methods, robust in-

Figure 2 :

 2 Figure 2: Declination, inclination, and intensity at Earth's surface from the International Geomagnetic Reference Field model (IGRF-12) for 2015 (Thébault et al., 2015) for three cases: (a) the axial dipole contribution alone; (b) the tilted dipole contribution alone; and (c) all dipole and non-dipole contributions.

Figure 3 :

 3 Figure 3: Declination map from IGRF-12 for 2015. White symbols mark Earth's rotation axis (open circle), magnetic pole (triangle), and geomagnetic pole (plus). Coloured symbols indicate three example sites (filled circles) and their respective VGPs (stars, see text for explanation) according to the present day field structure. The yellow star marks one of the lowest latitude VGPs in the present-day field. It is found in the southern Indian Ocean (62 • S, 62 • E), not seen on this projection. Blue is negative and orange is positive declination.

For a perfectly dipolar

  field, declination and inclination from globally distributed sites would result in the same VGP latitude and longitude, corresponding to the actual geomagnetic pole. Similarly, VDM values from intensities all over the globe would give the true dipole moment; this is only true for VADM values if the pure dipole field is not tilted away from the rotation axis. The influence of the non-dipole field results in differences and differing variations in VGP, VDM, and VADM depending on the location. In Fig. 3 we demonstrate three examples of how VGPs calculated from directional data at different sites deviate from the actual geomagnetic pole for the present day field (underlain is the declination map from IGRF 2015, as shown in Fig. 2c). For chronological applications, the concepts of VDM, VADM, and VGP are useful under the assumption that non-dipole variations are the same at the compared locations. More generally, they can be used to investigate the varying influence of non-dipole fields at different locations. Determining VADMs, VDMs, and VGPs from large numbers of palaeomagnetic results provides information about how well the GAD approximation holds. Large scatter in such data sets indicates strong non-dipole PSV, and systematic deviations of the mean VGP from the rotation axis point to persistent nonaxial-dipole contributions.

3 .

 3 Palaeomagnetic data and work flow Palaeomagnetic data come from two main sources: (i) thermal remanent magnetizations (TRMs) preserved in burnt archaeological materials (e.g., pottery) and igneous rocks which provide spot readings of the field in time and (ii) depositional remanent magnetizations (DRMs) acquired in sediments. Palaeomagnetic data obtained from non-clastic sedimentary rocks, such as flowstones and speleothems, which carry a chemical remanent magnetization, are a minor contributor to the global data set. The discrete or continuous nature of these recording mechanisms and the spatial and temporal distribution of archaeomagnetic, volcanic, and sediment data results in different insights into the palaeomagnetic field and, therefore, their application to palaeomagnetic dating.

  composed of stacked core data are also accommodated. An important aim of the GEOMAGIA50 sediment database is to allow sediment records to be re-evaluated using individual core data. This allows any researcher to remake the age-depth model for a core or composite and transfer this to the palaeomagnetic data. As of July 2018, GEOMAGIA50 contains ∼5600 Holocene archaeological and volcanic directional data and ∼4570 palaeointensities. The Holocene sediment data consist of ∼20,000 entries on core depth and ∼3000 stacked data. Both sediment and archaeomagnetic sections of the database are periodically updated when new data are published.The HISTMAG database[START_REF] Arneitz | The HISTMAG database: combining historical, archaeomagnetic and volcanic data[END_REF] combines historical observations of the past ≈500 years with archaeomagnetic and volcanic data of the past 50 kyr. It was developed based on two existing databases, which have both been updated with additionally compiled data. The historical data mostly come from the compilation by[START_REF] Jonkers | Four centuries of geomagnetic data from historical records[END_REF]. This database consisted of a set of files with search options provided by Fortran codes and was only available from the authors. It has been updated with 4160 historical data (of one to three field components) in HISTMAG. The archaeomagnetic and volcanic data in HISTMAG mostly come from the GEOMAGIA50.v3 database, updated with 183 additional records from Austria, Germany, and Poland, that were not included in GEOMAGIA50 at the time.URLs for all of these databases are given in Appendix B. It is of the utmost importance for future progress in chronological applications of archaeomagnetic and palaeomagnetic dating that researchers continue to upload details of their work into these databases, to provide a direct means for others to access their data and to assess their accuracy and significance in the context of new and evolving studies.

Figure 4 :

 4 Figure 4: Historical (a,d), archaeomagnetic (b,e) and palaeomagnetic sediment (c,f) data distribution in space (a,b,c) and time (d,e,f). Brown (top in a,b,c) are declination, orange (middle in a,b,c) are inclination, and green (bottom in a,b,c) are intensity data. Density of points in a) and b) gives an idea of how the regions are covered by data in time. Sediment records in c) often sample the whole time interval (f). Distribution in time is given by both linear (top in d,e,f) and logarithmic (bottom in d,e,f) scale for better comparison of the vastly varying historical and archaeomagnetic/volcanic data. Note the different scale in the upper panel in d) compared to e) and f). Historical data are shown in 10-year bins in d), and archaeo-and palaeomagnetic data are given in 200-year bins in e) and f).

Figure 5 :

 5 Figure5: Overview of geomagnetic field models and dipole reconstructions for historical times through to 50 ka. Model references: gufm1[START_REF] Jackson | Four centuries of geomagnetic secular variation from historical records[END_REF], HHK1998[START_REF] Hongre | An analysis of the geomagnetic field over the past 2000 years[END_REF], SCHA.DI.00(Pavón-Carrasco et al., 2008b), SCHA.DI.00-F(Pavón- Carrasco et al., 2008a), GMADE2K.1[START_REF] Lodge | Towards a new approach to archaeomagnetic dating in Europe using geomagnetic field modelling[END_REF], VHLP2008[START_REF] Valet | Secular variation of the geomagnetic dipole during the past 2000 years[END_REF], CJL2000[START_REF] Constable | Global geomagnetic field models for the past 3000 years: transient or permanent flux lobes? Phil[END_REF], A FM, ASD FM, and ASDI FM[START_REF] Licht | Ensembles of low degree archeomagnetic field models for the past three millennia[END_REF], AmR[START_REF] Sanchez | Modelling the archaeomagnetic field under spatial constraints from dynamo simulations: a resolution analysis[END_REF], COV-ARCH and COV-LAKE[START_REF] Hellio | Time-correlation-based regression of the geomagnetic field from archeological and sediment records[END_REF], CALS3k.x, ARCH3k.1, and SED3k.1 (Korte and Constable, 2003, 2005; Korte et al., 2009; Korte and Constable, 2011), SCHA.DIF.3k (Pavón-Carrasco et al., 2009), CALS7k.2 (Korte and Constable, 2005), ARCH-UK.1 (Batt et al., 2017), ArcheoInt VADM (Genevey et al., 2008), SCHA.DIF.8k (Pavón-Carrasco et al., 2010), DE FNB, DE FNBK, and DE FNBKE[START_REF] Nilsson | Holocene geocentric dipole tilt model constrained by sedimentary paleomagnetic data[END_REF][START_REF] Dunlop | Physical basis of the Thellier-Thellier and related paleointensity methods[END_REF]) pfm9k.1, pfm9k.1a, and pfm9k.1b (Nilsson et al., 2014), OH1992[START_REF] Ohno | Geomagnetic poles over the past 10,000 years[END_REF], OH1993(Ohno and Hamano, 1993), CALS10k.x (Korte et al., 2011;[START_REF] Constable | Persistent high paleosecular variation activity in Southern Hemisphere for at least 10 000 years[END_REF], HFMx, and HFM.OL1.A1[START_REF] Panovska | Limitations in paleomagnetic data and modelling techniques and their impact on Holocene geomagnetic field models[END_REF][START_REF] Constable | Persistent high paleosecular variation activity in Southern Hemisphere for at least 10 000 years[END_REF], YOS2000[START_REF] Yang | Variations in the geomagnetic dipole moment over the last 12000 years[END_REF], SHA.DIF.14k(Pavón-Carrasco et al., 2014b), MS1982 (McElhinny and Senanayake, 1982), and GEO-MAGIA VADM[START_REF] Knudsen | Variations in the geomagnetic dipole moment during the Holocene and the past 50 kyr[END_REF].

Figure 6 :

 6 Figure 6: Predictions from seven models for regions with different data coverage (see text and Fig. 4 for details): a) central Europe (50 • N, 15 • E), b) east Asia (36 • N, 128 • E), c) the southern tip of South America (55 • S, 68 • W), and d) South Africa (31 • S, 23 • E). Solid lines are from models with both archaeo-and palaeomagnetic data, including sediment records (see model names in panel b). Dashed lines are from models based only on archaeomagnetic and volcanic data. Model uncertainties are not shown (to make the figure clearer).

  cordant. For example, in the cases shown in Fig.6the model uncertainties are generally on the orders of up to 7 • in declination, 6 • in inclination and 6 µT in intensity. Larger uncertainties are seen for Southern South America in some models, but uncertainties for South Africa are not systematically larger than for Europe in any of the models. Therefore it is important to consider how well-constrained a global model is for a given region and time interval when using it for palaeomagnetic chronology refinement. Two recent approaches aim at better characterizing model uncertainties.[START_REF] Sanchez | Modelling the archaeomagnetic field under spatial constraints from dynamo simulations: a resolution analysis[END_REF] utilized prior information from mean and covariance from an ensemble of magnetic field states of a geodynamo simulation and presented a series of global snapshot models as a step toward archaeomagnetic data assimilation into physics-based models.[START_REF] Hellio | Time-correlation-based regression of the geomagnetic field from archeological and sediment records[END_REF] proposed a time-correlation based regression using a priori information based on spatial and temporal characteristics of the present-day field. They concluded that their models have statistically coherent uncertainties, and that the model including lake sediments is much better constrained (due to the much larger number of data) than the the one based only on archaeomagnetic data. However, their model predictions, like other global models, do not fully resolve centennial field variations. Many models are freely available and the URLs are given in Appendix B.

  important role in failure to refine individual ages, particularly if only one or two field components have been recovered. Moreover, variations in reference curves and models might be significantly smoothed compared with the amplitudes recovered by archaeomagnetic or volcanic spot values. Sediment records in principle should reflect the same field variations as a reference record and, thus, a direct correlation should be straightforward. However, depending on sedimentation rate and the chosen reference curve, either the record or the curve might be notably smoother than the other, which complicates direct correlation. Moreover, the magnetization lock-in age obtained from such a correlation might differ from the sediment age. Uncertainties both in recovered magnetic field values and in reference curves should not be neglected, and it should be considered that estimates of both model and data uncertainties might underestimate true errors.

  developed the REN-DATE software to perform these tasks and Pavón-Carrasco et al. (2011) subsequently implemented it in a freely available, easy-to-use Matlab routine (see URL in Appendix B). An example of its application is shown in Fig. 7, where sample values are the average results for pottery kiln Serdica 1 investigated by Kovacheva et al. (2004) and the dating curves come from the SHA.DIF.14k model. The resulting most probable age interval ranges from 1599 to 1680 AD and overlaps well with the original result obtained with REN-DATE from the Bulgarian reference curve of 1513 -1683

Figure 7 :

 7 Figure 7: Example of archaeomagnetic age refinement using the Matlab tool of Pavón-Carrasco et al. (2011). Output model curves for declination, inclination, and intensity with uncertainty estimates (top, red lines) with the field value of the material to be dated as blue line with green uncertainty limits. The middle row of panels are the individual component probability density functions with probable ages above the 95% significance level (green line) shaded in grey. Bottom panels are combined probability density function results next to the location map.

  two cases they found good agreement with the archaeological expectation, but in four cases ages were younger than expected and two examples could not be dated. Similarly,[START_REF] Morales | Archeomagnetic dating of some Pre-Columbian pottery fragments from Northern Mesoamerica: Implications for the chronology of Central Mexico during the Epiclassic period[END_REF] used intensity variations from a 3 kyr model to determine the ages of Mexican pottery and found one result to be in excellent agreement with archaeological expectation, while another was lower than any part of the curve. It remains unclear whether the failed dating attempts in both examples are due to limited resolution of the curves, earlier ages of the artefacts, or inaccuracies in archaeomagnetic results.The effect of insufficient reference curve resolution, which in particular can be a problem for curves obtained from global models, is shown by[START_REF] Stillinger | Archaeomagnetism as a complementary dating technique to address the Iron Age chronology debate in the Levant[END_REF], who developed a regional Near East reference curve that better describes a strong intensity maximum around 1000 BCE than any existing global model. In a study of material from hearths in South Korea,[START_REF] Shin | On the validity of archeomagnetic dating method in Korea: a case study[END_REF] found good overlap with radiocarbon ages when using a regional Japanese curve and a 3 kyr global model curve for archaeomagnetic dating of directional data. There was no systematic relationship between the ages obtained from the archaeomagnetic and radiocarbon methods, so they concluded that neither of the two PSV curves accurately represents past magnetic field variations in Korea.Several examples of successful age refinement of volcanics by palaeomagnetic dating are reported in the recent literature. For instance,[START_REF] Speranza | Holocene eruptive history of the Stromboli volcano: constraints from paleomagnetic dating[END_REF][START_REF] Speranza | Paleomagnetic dating of the most recent silicic eruptive activity at Pantelleria (Strait of Sicily)[END_REF];[START_REF] Tanguy | New archeomagnetic and 226 Ra-230 Th dating of recent lavas for the Geological map of Etna volcano[END_REF] concluded that palaeomagnetic dating of directional results from Stromboli, Pantelleria, and Etna lava flows in Italy, respectively, considerably narrowed ages obtained from radiometric methods over the whole Holocene using regional reference curves. Similarly, the eruptive histories of two Mexican volcanic areas were refined by[START_REF] Böhnel | Palaeomagnetic dating of two recent lava flows from Ceboruco volcano, western Mexico[END_REF];[START_REF] Mahgoub | Paleomagnetically inferred ages of a cluster of Holocene monogenetic eruptions in the Tacámbaro-Puruarán area (Michoacán, México): implications for volcanic hazards[END_REF], who used full vector data and reference curves from a global field model. Most of these studies also point out how chronological information from different sources can complement each other to determine the most probable ages. However, the study byRoperch et al.

(

  2010) used a global model to refine the chronology of a record from the Beaufort Sea constrained by only one radiocarbon age. Using only declination for correlation to the model curve and by comparing directional results to other Arctic and North American records, they concluded that millennial-scale directional PSV features recovered by the global model in the region allowed palaeomagnetic age refinement. Suteerasak et al. (2017) successfully dated three Bothnian Bay (between Sweden and Finland) cores through compari-

  ferences, we show examples of VADM, declination, and inclination anomaly from the CALS10k.2 field model at increasing distances from a central, midlatitude location (Windermere, UK) in Figs. 8, 9, and 10. It is obvious how

Figure 8 :

 8 Figure 8: CALS10k.2 model predictions for VADM (intensity) at Windermere, United Kingdom (54 • N, -4 • W) (black lines), plotted against model predictions to the N, S, E, and W of the signal at 500, 1000, 2000, and 4000 km distances (coloured dots on projections to the right indicate the location of the time series).

Figure 9 :

 9 Figure 9: CALS10k.2 model predictions for declination ( • ) at Windermere, United Kingdom (54 • N, -4 • W) (black lines), plotted against model predictions to the N, S, E, and W of the signal at 500, 1000, 2000 and 4000 km distances (coloured dots on projections to the right indicate the location of the time series).

Figure 10 :

 10 Figure 10: CALS10k.2 model predictions for ∆I ( • ) at Windermere, United Kingdom (54 • N, -4 • W) (black lines), plotted against model predictions to the N, S, E, and W of the signal at 500, 1000, 2000 and 4000 km distances (coloured dots on globe plots to the right indicate the location of the time series).

  . For sedimentary records, which are inherently smoothed due to the continuous nature of sedimentation, comparisons of overall dipole-dominated field model predictions are likely also valid, and correlation distances shown in Figs. 8, 9, and 10 might be applicable depending on the resolution of the record. Considerable discrepancy might be expected, on the other hand, from point-value field measurements and smoothed model predictions. Comparisons with surface field model predic-

Figure 11 :

 11 Figure 11: Maps of the correlation coefficients of declination, inclination, and intensity time series from CALS10k.2 compared with the CALS10k.2 prediction for Windermere (yellow star). Contour lines for R=0.8 (dashed yellow) and R=0.9 (solid white) are shown for reference.

Figure 12 :

 12 Figure 12: Temporal evolution of the dipole latitude from numerical dynamo simulation Case A. The time axis scales with magnetic diffusion times, where 1 MDT is 40,868 calendar years. Green bars mark a 45• dipole tilt, and give a threshold for field excursions to be considered global. The red-shaded interval indicates the time span used to provide synthetic data for the correlation analysis.

Figure 13 :

 13 Figure13: Maps of R95, the 95% confidence limit of R, for inclination, declination, intensity, and VGP latitude time series from Case A dynamo simulation compared to a reference time series also from Case A at latitudes 0, 45, and 90 • N/S (locations shown by yellow star). Contour lines for R95 = 0.8 (dashed yellow) and R95 = 0.9 (solid white) are shown for reference. Also shown, on the right-hand panel, are plots of the average distance in km from the reference point to the contour lines for R95 = 0.8 and R95 = 0.9 shown in the corresponding maps.

  VGP path moves clockwise, the inclination maximum at site 2 occurs after the maximum for site 1. Depending on the rate of VGP change, this results in some amount of temporal offset. The difference between inclination and declination at different sites, although the VGP has the same coordinates, explains why VGP variations in Fig.13are correlated over a greater distance than inclination and declination alone. It also highlights that even with a dipolar configuration, offsets in the timing of inclination and declination are example of the relationship between VGP variations and declination and inclination changes assuming a dipole field in which a VGP with constant latitude varies between 180 • and -180 • longitude. (a) VGP path (rainbow line) at a VGP latitude of 75 • . Site 1 (grey circle) and site 2 (black circle) are two locations 45 • apart for which declination and inclination changes are shown in (b) and (c), respectively. The rainbow colour scales above (b) and (c) correspond to the rainbow band in (a).

  VGP positions are shown for a portion of the MD99-2269 record of Stoner et al. (2007). In Fig. 15b VGP latitude change with time is shown for the same record. To demonstrate how changing the observation site influences declination and inclination when assuming that the observed directional variations originate from a tilted geocentric dipole field, we use VGPs from MD99-2269 to calculate the expected declination and inclination for four locations at increasing distances from MD99-2269 (Fig. 15c, d). At Hvítárvatn, a lake in central Iceland used in the PSV study of Ólafsdóttir et al. (2013) (approximately 200 km or 1.8 • from MD99-2269), the relocated inclination and declination are reproduced almost exactly, with insignificant magnitude differences for the directions. At 30 • (∼1340 km) due west of MD99-2269, recalculated directions are broadly similar to those at the original site. Some directional variations can clearly be correlated; however, the magnitude and duration of others are now quite different, e.g., between 2.5 and 3 ka, and 5 and 6 ka. At 45 • (∼2010 km) due west of MD99-2269 correlation to the original record is ambiguous for most of the time series. Some features are similar in either inclination or declination, but offset by around 200 years. At 90 • (∼4030 km) due west of MD99-2269, there is little resemblance between the original and recalculated records. These differences stem purely from how a tilted dipolar field would be observed at different observation sites. For example, between 2.5 ka and 3 ka, the poles are at their most southerly latitudes (over central Europe) and are far from the site 90 • west of MD99-2269. The dipole is therefore tilted away from the observation site and the corresponding inclinations are shallow. Conversely, at MD99-2269, VGPs are much closer to the site and inclinations are, therefore, steeper.

Figure 15 :

 15 Figure 15: Illustration of the dependence of declination and inclination on the observation site for a dipolar field. (a) Small coloured circles: VGPs calculated from the 2-6 ka portion of the directional record of MD99-2269 (Stoner et al., 2007). The different colours denote 0.5 kyr segments of the record (see (b) for the time divisions). The large black circle is the site of the MD99-2269 record; the large grey circle is Hvítárvatn, Iceland; the large coloured circles are observation sites at certain degrees longitude away from MD99-2269. (b) VGP latitude from the same portion of the MD99-2269 shown in (a); colours correspond to those in (a). (c) Declination and (d) inclination at the observation sites shown in (a). The data shown for Hvítárvatn, 30 • , 45 • and 90 • west of MD99-2269 were calculated (relocated) using the declination and inclination data of MD99-2269 assuming the observed variation resulted purely from variations of a geocentric, tilted dipole field. Coloured bands above the plots correspond to the VGP circles in (a) and the VGP latitude segments in (b).

  6.4. Temporal correlation uncertainties Temporal errors resulting from PSV correlations are difficult to quantify because no standardized methods for PSV correlation or tie point selection have been determined. Comparing model curves from different locations by wiggle-matching (i.e., one-to-one matching of specific patterns) can introduce large age errors, as illustrated below with two examples. To illustrate the introduction of age errors by wiggle matching, we use CALS10k.2 model predictions for different locations and field components with a correlation coefficient of R ≈0.8 with the Windermere model curve, (compare Fig. 11) and we compare them with that curve. The records are matched in Fig. 16 via dynamic programming (the Match-2.2 algorithm of

Figure 16 :

 16 Figure 16: Results of matched VADM (top), inclination anomaly (middle), and declination (bottom) time series from the CALS10k.2 model for the mid-latitude Windermere location with targets (different for each component) from locations that show a correlation of R≈0.8. Left-hand panels are normalized field components before running the Match-2.2 algorithm of Lisiecki and Lisiecki (2002). Right-hand panels are results of the matching. Mean offsets in years for the time series are given in the right-hand panels.

Figure 17 :

 17 Figure17: Illustration of correlation from numerical dynamo simulation runs for 6,000 years (based on magnetic diffusion time) during non-excursional field behavior. Signals and targets were chosen for a correlation coefficient of R≈0.8. Left-hand panels are normalized field components intensity, inclination, and declination for signals and targets before running the Match-2.2 algorithm of[START_REF] Lisiecki | Application of dynamic programming to the correlation of paleoclimate records[END_REF]. Centre panels are results of the match, and right-hand panels are offsets in years at each time step. Note that components have been matched individually.

  . Individual uncertainties are higher than previously considered, with a wide range of values indicating that each record and component is different in quality. Median values determined for 73 Holocene sediment records are 5.9 • for inclination, 13.4 • for declination, and 11 µT for the standardized RPI based on calibration with the CASL7K.2 field model. For each RPI record, the uncertainty must be converted to absolute values with the same calibration factor as for the RPI (see Sec. 3.4.2).

	particular, rather than supplying data targeted for specific kinds of analy-
	ses it provides a searchable compilation and archive of published data and
	allows users to select data by criteria based on wide-ranging metadata. It
	is the user's responsibility to determine, on the basis of published work and
	archived metadata, whether the data may be suitable for any given purpose.
	Given sufficient fundamental measurement information it can be possible to
	produce new interpretations at levels ranging from the laboratory to global
	modelling. Tauxe et al. (2016) provide a software package to facilitate (re-
	)interpretation of laboratory experiments and show how specific work-flows
	and interpretation tools can be used to generate data in formats suitable
	for direct upload of palaeomagnetic and rock magnetic results into MagIC.
	The MagIC (Magnetics Information Consortium) database is a compre-
	hensive and flexible sample-based database designed to accommodate all
	kinds of palaeomagnetic, and rock magnetic and related data and meta-
	data such as location, materials, and chronological information. The data
	range from basic summaries of published results to individual experimental
	measurements. MagIC differs from GEOMAGIA and HISTMAG (described
	below) in a number of ways. The structure of the data archive is linked to δ 18 O, tephra, varve thermoluminiscence (TL) and optically stimulated lu-
	data produced and/or analysed in individual (usually peer-reviewed) pub-miniscence (OSL) ages) to be entered at different analysis levels. Data from
	lications. It is not restricted to a specific time interval, but spans the age individual sediment cores can be entered for core depths or composite depths,
	of the Earth. The design of MagIC is driven by community need to ad-with associated inferred age from the published age-depth model. Records
	dress broad scientific challenges through standardizing diverse datasets and
	maintaining an open archive for published rock and palaeomagnetic data. In

3.6. Data reporting

Published palaeomagnetic data for dating can best be used if they are made available through publicly accessible online databases. The need for global data compilations was recognized several decades ago, and the International Association for Geomagnetism and Aeronomy (IAGA) encouraged the development of a series of palaeomagnetic databases from the late 1980s on (see, e.g.,

[START_REF] Constable | Centennial-to millennial-scale geomagnetic field variations[END_REF]

. The past 15 years has seen the evolution of previous efforts into several new data initiatives, each developed with specific goals in mind. The current versions, which are relevant for Holocene times, are MagIC, GEOMAGIA, and HISTMAG. All three have interactive online search interfaces that allow data selection by various criteria ranging from region or age to details such as material, dating method, and many more given by the metadata. Flexibility in the specific details also allow individual researchers to develop their own strategies for uploading data. As of August 2018, 49 researchers have made 4,295 contributions to MagIC (each associated with a publication). A long term goal is for individuals to upload their own data as they are published, but the first upload to MagIC has often not been initiated by the publication's author(s). The level of detail supplied in various contributions spans a broad range. Over 200,000 sites are represented worldwide and more than 4.6 million individual measurements. MagIC has substantial overlap with GEOMAGIA50, although there are distinct records that are unique to each database.

GEOMAGIA50 is an online database of palaeomagnetic and chronological data from archaeological materials, volcanic rocks, and sediments spanning the past 50 kyr. Details of the database are give by

Brown et al. (2015b)

and

Brown et al. (2015a)

. Data from archaeological materials and volcanic materials are reported at the site level (e.g., archaeological horizon or lava flow).

Data are organized by age and the data reported include palaeomagnetic directions, α 95 , Fisherian concentration parameter, k, palaeointensity and uncertainty, metadata on the geographic location of sampling sites, palaeomagnetic methods (direction and intensity), and the types of specimens measured. The sediment database allows palaeomagnetic (directional and RPI data), rock magnetic (e.g., ARM, IRM, k), and geochronological data ( 14 C,

  used for age-depth model construction should be published with the magnetic data, including but not limited to: core and composite depth of age tie-points, depths of uncalibrated radiocarbon ages with 1σ uncertainties, calibrated radiocarbon ages with 2σ uncertainties, calibration/reference curves used, software used to create age-depth models, interpolation method, rejected or outlying ages, reservoir effects and associated errors, or additional age information used (e.g., OSL dates or tephra ages are sometimes used as additional tie points in 14 C chronologies). When an age-depth model is refined by PSV correlations, tie points for these correlations and the reference curve used to identify tie points should be provided such that the age model can be reconstructed independently without PSV tie points. This ensures age-depth model reproducibility and usefulness of the data as inputs for improved reference curves and models.

	Many open-source palaeomagnetic databases like MagIC and GEOMA-
	GIA50.v3 (Appendix B) accommodate palaeomagnetic, archaeomagnetic,

  are partly caused by persistent structures in the field and are specific to how the field varied at and around Windermere over the past 10,000 years. The spatial correlation maps and the examples in Figs. 8 to 10 suggest that intensity records are most suitable for age refinements over large distances. This is due to

	the strong influence of (global) dipole strength on intensity variations. The
	high-frequency nature of PSV features in inclination and declination com-
	ponents, on the other hand, indicate that correlation of these components
	on shorter timescales or at shorter correlation distances may produce higher
	precision results (Figs. 9, 10). However, due to the strong smoothing in
	the model, all predictions likely over-estimate correlation distances, and cor-
	relation distances vary greatly with location and with the field component
	analysed.
	Moreover, correlation between real datasets and surface model predictions
	may introduce age offsets in records where local non-dipole field features are
	detected but are matched to dipolar features in surface models. Model predic-
	tions from areas of lower data density (e.g., high latitudes and the Southern
	Hemisphere in CALS10k.2) may only have longer wavelength field structures,

dependent ages or chronologies, and diligent assessment of uncertainties for both palaeomagnetic data and ages. Palaeomagnetic databases that include metadata on methods and all relevant information to update ages if additional information becomes available provide a valuable basis for improving reference curves and global models and we strongly encourage the submission of any new results to such databases.

Acknowledgements and author contributions

Ideas and the outline for this review paper were developed collaboratively by all authors at the 8 th Nordic Palaeomagnetism Workshop held in Leirubakki, Iceland, from 30 th September to 7 th October 2017. Correlation analyses were carried out by AN and SG, dynamo simulation data were provided by IW, and MK coordinated the writing of the manuscript, to which all authors contributed. All results and presented material were discussed among all authors. We thank Joseph Stoner for sharing the MD99-2269 sediment record. Figures 2,3, and 4 were produced with free software Generic Mapping Tools [START_REF] Wessel | New, improved version of Generic Mapping Tools released[END_REF]. Figure 6 was produced with free software plotxy by Bob Parker (https://igppweb.ucsd.edu/ parker/Software/). Andrew Roberts and an anonymous reviewer a thanked for constructive comments on the initial mansucript.

Appendix A. Geodynamo simulation

In this study, we numerically simulated temporal geodynamo behaviour over ≈ 2.4 Myr, with a temporal resolution of about 60 years. The geodynamo simulation used here, case A, represents a self-sustaining, convecting magnetic dynamo and is simulated using the MagIC code [START_REF] Wicht | Inner-core conductivity in numerical dynamo simulations[END_REF].

The code provides a full three-dimensional numerical solution of a set of largely differ from theoretical values for Earth, temporal behaviour from the simulation has Earth-like variations, such as reversals and excursions paired with periods of stable field polarity (chrons). The numerical simulation shows a reversal rate which is comparable to the observed reversal rate of the geomagnetic field. Here, we consider a polarity change to be a field reversal if the periods of oppositely oriented stable field before and after lasted at least one magnetic diffusion time (≈ 40 kyr). In its total run time of ∼ 2.4 Myr, the numerical simulation has at least 7 such reversals, which is similar to the reported geomagnetic field reversal rate for the last 5 Myr [START_REF] Cande | Revised calibration of the geomagnetic timescale for the late Cretaceous and Cenozoic[END_REF][START_REF] Constable | On the rates of geomagnetic reversals[END_REF]. The field morphology is also broadly similar to the geomagnetic field. [START_REF] Christensen | Conditions for Earth-like geodynamo models[END_REF] suggested various criteria that characterize field morphology in the dynamo simulation. A primary property of the geomagnetic field is the dominance of the axial dipole: its dipolarity.

For Earth, the dipolarity is d = 1.4 [START_REF] Christensen | Conditions for Earth-like geodynamo models[END_REF], the value for our simulation is d = 0.4, providing more small-scale field variations as desired here. Another criterion rates the equatorial symmetry of the field, this is s = 1.0 for Earth, and a very similar value of s = 0.97 for our simulation.

Both values indicate a weak preference for equatorial anti-symmetry. Magnetic field zonality, also used to characterize dynamo simulations, is lower with z = 0.03 in case A compared to z = 0.15 for Earth.

Appendix B. Digital resources

In the following we list the URLs of freely accessible databases, computer software and global geomagnetic field models that are mentioned in this review paper and can be useful for archaeo-and palaeomagnetic age refinement.

Appendix B.1. Databases and archives useful to Holocene palaeomagnetic age refinement

• GEOMAGIA50.v3: http://geomagia.gfz-potsdam.de/index.php

Database containing archaeomagnetic and volcanic data and palaeomagnetic records for the last 50 kyr with detailed metadata, particularly including dating information (Brown et al., 2015a,b). Online query forms exist for archaeomagnetic/volcanic and for sediment data.

• MagIC: https://www.earthref.org/MagIC Digital data archive for rock and palaeomagnetic data with portals that allow users access to archive, search, visualize, and download data.

• HISTMAG: http://www.conrad-observatory.at/zamg/index.php/dataen/histmag-database Database combining historical, archaeomagnetic, and volcanic data [START_REF] Arneitz | The HISTMAG database: combining historical, archaeomagnetic and volcanic data[END_REF]. Registration is required to use the web interface.

• PANGAEA: https://www.pangaea.de/ Broader data archive that includes sediment records.

Appendix B.2. Dating tools and global model PSV curves

• Matlab dating tool: https://earthref.org/ERDA/1134/ Matlab code for determining age probability densities for input of one to three geomagnetic field components from a variety of models and reference curves [START_REF] Pavón-Carrasco | A Matlab tool for archeomagnetic dating[END_REF].

• REN-DATE archaeomagnetic dating software: http://dourbes.meteo.be/aarch.net/onlytxt/rendate.otxt en.html Program underlying the above Matlab code for Microsoft Windows operating systems [START_REF] Lanos | Bayesian inference of calibration curves: application to archaeomagnetism[END_REF]).

• GEOMAGIA50.v3: http://geomagia.gfz-potsdam.de/index.php Database query forms allow calculation of PSV curves from several global models for any location on Earth (Brown et al., 2015a,b). Direct link to download these models [START_REF] Licht | Ensembles of low degree archeomagnetic field models for the past three millennia[END_REF] and Matlab code to obtain PSV curves.

Appendix B.3. Palaeomagnetic and correlation utility programs

• Match: http://www.lorraine-lisiecki.com/match.html Software package using dynamic programming to find the optimal alignment of two signals using penalty functions to constrain sediment accumulation rates, available as a command line version or with a Matlab interface [START_REF] Lisiecki | Application of dynamic programming to the correlation of paleoclimate records[END_REF].

• PmagPy: https://earthref.org/PmagPy Software package developed for palaeomagnetic data analysis written in Python [START_REF] Tauxe | PmagPy: Software package for paleomagnetic data analysis and a bridge to the Magnetics Information Consortium (MagIC) Database[END_REF]