Calcite Twin Formation, Measurement and Use as Stress–Strain Indicators: A Review of Progress over the Last Decade

: Mechanical twins are common microstructures in deformed calcite. Calcite twins have been used for a long time as indicators of stress/strain orientations and magnitudes. Developments during the last decade point toward signiﬁcant improvements of existing techniques as well as new applications of calcite twin analysis in tectonic studies. This review summarises the recent progress in the understanding of twin formation, including nucleation and growth of twins, and discusses the concept of CRSS and its dependence on several factors such as strain, temperature and grain size. Classical and recent calcite twin measurement techniques are also presented and their pros and cons are discussed. The newly proposed inversion techniques allowing for the use of calcite twins as indicators of orientations and/or magnitudes of stress and strain are summarized. Beneﬁts for tectonic studies are illustrated through the presentation of several applications, from the scale of the individual tectonic structure to the continental scale. The classical use of calcite twin morphology (e.g., thickness) as a straightforward geothermometer is critically discussed in the light of recent observations that thick twins do not always reﬂect deformation temperature above 170–200 ◦ C. This review also presents how the age of twinning events in natural rocks can be constrained while individual twins cannot be dated yet. Finally, the review addresses the recent technical and conceptual progress in calcite twinning paleopiezometry, together with the promising combination of this paleopiezometer with mechanical analysis of fractures or stylolite roughness.


Introduction
A twin is a planar crystallographic defect that separates two homogeneous parts of a single crystalline species in a symmetric way [1] (Figure 1A-C). In the classical theory of deformation twinning, the original (parent) lattice is re-orientated by atom displacements which are geometrically equivalent to a simple shear of the lattice ( Figure 1C).
At temperatures below those at which individual atoms are mobile, slip and twinning are the dominant deformation mechanisms that enable a solid to change its shape under the action of applied stress. In calcite, the most commonly observed deformation twinning rule corresponds to twinning along {01-12} (or {018}) e-planes. Deformed calcite readily acquires intracrystalline twin deformation lamellae, caused by lattice rotation of 38.2 • between rows of atoms parallel to the e-plane in response to shear stress [2]. Due to the trigonal symmetry of calcite, three symmetrically equivalent e-twins may occur ( Figure 1C).  . The twin planes e1, e2 and e3 are shown as great circles. Each twin plane contains the twinning direction shown as a point. For each twin plane, the arrow is parallel to the twinning direction; its head indicates that the upper part of the crystal moves upward, toward the C-axis, as a reverse microfault. (D) Geometry of calcite twinning. C and C' refer to the optical axes of the host and twinned portion of the crystal, respectively. e is the pole to the twin plane and g is the twinning direction. (E). Photographs of twinned calcite grains in a vein in natural light (left) and under cathodoluminescence (right). Note that twin lamellae sometimes display a distinctive signature with respect to host grain.
Low-temperature calcite plasticity (0-300 • C) corresponds to the dominance of etwinning, making e-twins the most common microstructures in limestones ( Figures 1A,B  and 2A,B). Once a twin is formed, progressive deformation can create new twins elsewhere or deform the twinned rock portion by other mechanisms. Mechanical twinning can however accommodate only a fraction of the bulk strain affecting the rock. Consequently, twinning deformation in upper crustal carbonate rocks generally co-exists with other  Calcite twins have given rise over years to two main kinds of investigations. On the one hand, calcite twinning has been studied as one crystal plastic deformation mechanism through deformation experiments to better constrain factors that control twin formation and the contribution of twinning to the deformation of calcite grains at various conditions. On the other hand, calcite twins have been widely used as a paleogeothermometer [4,5] as well as strain [6,7] and stress (e.g., [3,8]) gauges in tectonic studies. may be either to produce the very high stresses required for homogeneous nucleation or to form nuclei more directly by re-arrangements in their core structures [11]. Elastic twins nucleated at such stress heterogeneity on a grain boundary propagate across the grain, probably very rapidly with increasing stress. It is likely that because of neighbouring grain constraints, the number of twin planes with positive Schmid factors available in a grain influences the easiness of this elastic process (Covey-Crump, pers. comm., 2021). Observations in calcite single crystals [22] and neutron diffraction experiments in polycrystalline limestone [17,18] show that elastic twins may nucleate and grow even when the applied stress is below the aggregate yield point, i.e., deformation is reversible.
At loads greater than the aggregate yield point, inelastic (plastic) permanent twinning occurs, and twin lamellae thicken owing to repulsive interactions between twin dislocations and matrix defects, e.g., cracks, dislocations, and other twins, all of which may hinder the reverse motion of twin boundaries [23]. It is possible, but still unclear if elastic twinning can be really seen as some kind of early stage of development of a permanent twin, i.e., some precursor twin.
In the experiments by Parlangeau et al. [13], the produced macroscopic stress-strain curves exhibit three major stages ( Figure 3A). The first two stages I and II correspond to a preloading stage and a pseudo-elastic stage, respectively. The elastic-plastic transition is then reached and stage III corresponds to the onset of plastic deformation twinning. The plastic yield strength corresponds to twinning activation and is evidenced by the change in slope of the stress-strain curve ( Figure 3A). Beyond this threshold, deformation stage IIIa is characterized by strain hardening, which is related to the early onset of crystal plasticity with the activation of the very first isolated mechanical twins, followed by slight densification of twin lamellae until 1-2% shortening. In stage IIIb, twinning-related plastic deformation occurs with a strong increase in twin density and bursts of thickening of twin lamellae. The corresponding macroscopic stress-strain curve is strongly serrated, reflecting episodically fluctuating crystal strength, with strain hardening increments alternating with sudden stress drops. Sometimes, the stress-strain curve exhibits a low-stress plateau indicating nearly steady-state flow stress ( Figure 3A). Most irreversible strain is achieved during stage IIIb. Interestingly, twinning-dominated deformation is accompanied throughout the experiments by micro-cracking that allows for the accommodation of local strain incompatibilities related to anisotropic twinning deformation. The observation that micro-cracking and twinning in single crystals probably go hand in hand can be extended to calcite polycrystals deformed at room pressure [18].
A major contribution of this study was to establish that from a phenomenological point of view, thickening of twin lamellae at low temperatures occurs in a sequence comprising (1) densification of twins, (2) thickening of individual lamellae and (3) merging of thickened lamellae ( Figure 3B). The widest twin lamellae appear to be the final result of several merging events of thinner twin lamellae [13]. The growth of twins by the nucleation of twinning dislocations on planes parallel and contiguous to the coherent twin boundary provides a viable mechanism for twin widening [11]. The experimental condition did not enable [13] to reach the subsequent stage of increased hardening where twinning has saturated and other crystalline deformation mechanisms must activate.
The occurrence of intermittent stress drops during deformation of single crystals beyond the yield point ( Figure 3A) is correlated with twinning events, i.e., the appearance of new twins or twin thickening [13,24] (Figure 3), which suggests that the stress released by, and possibly associated with, gliding of twin dislocations related to twin propagation and widening is small compared with that necessary to 'nucleate' them. This is in line with the idea that it may be necessary to consider separately the stress needed for nucleation of a twin and the (usually lower) growth stress [25].   [13] and showing the three stages discussed in the text. (B) High-resolution pictures showing twin development (stage III) from the formation of thin twins to late twin thickening. Each photograph has a size of 3 × 6 mm. The twin lamellae are represented in light green, and a high density of twin lamellae is highlighted using dark green. Fractures are displayed in pink. See text and [13] for details.

The Concept of CRSS for Twinning
One of the pending questions pertaining to e-twinning is related to the existence of a critical resolved shear stress (CRSS). If it does exist, what are its actual significance, value and possible dependence on temperature, strain, strain rate and grain size? These questions are of primary importance for tectonic and paleopiezometric studies. Indeed, except a couple of paleopiezometric approaches [18,26], all methods that aim at evaluating differential stress magnitudes from randomly oriented grains in an aggregate assume that mechanical twinning occurs once the shear stress across one of the three e-twin planes and resolved along the twin vector reaches (overcomes) the CRSS for twinning. Such a law would be analogous to Schmid's law for slip.
The strength of a sliding system (e-twinning or r-or f-gliding) is conventionally expressed by a CRSS. The CRSS corresponds to the resolved shear stress along the sliding plane along the sliding direction that must be reached to induce a significant plastic (permanent) deformation, i.e., to induce motion of a number of dislocations, so that sliding becomes macroscopically observable independently of the orientation of the deformed grain. Such behaviour is commonly associated with a critical point on the stress-strain curve for a monocrystal. The value of the CRSS is given by: τ C = σ × S, where σ corresponds to the applied differential stress at the critical point and S is the Schmid's factor, such as S = cos α × cos β, with α the angle between compression and the normal to the twin plane and β the angle between compression and the twin vector. The resolved shear stress along the twin vector is maximum when α and β are equal to 45 • , S varying between −0.5 and 0.5 depending on crystal orientation.
Experimental calibration of the stress at which the twins actually formed is inherently difficult. A fundamental question that therefore arises when trying to understand the onset of twinning is whether stress is homogeneous at the grain scale so that the magnitudes (and orientations as well) of the grain scale stresses reflect those applied to the bulk sample, or, alternatively, if a grain twins as a result of the local stress state within the microstructure. In other words: does a twin initiate in response to local/internal or external stresses? Internal stress distribution within a rock is intrinsically inhomogeneous [27] and the relationship between the force acting on each calcite grain and the bulk applied force remains uncertain.
If the picture for nucleating and propagating a twin across a grain from some stress heterogeneity is correct, and that (a) the nature of the local heterogeneity (affecting the stress required to activate it), and (b) the ease with which the strain can be accommodated given neighbouring grain constraints, are important, the existence of a CRSS is suspect at first glance [18,28]. Two perspectives have been discussed by [18]. The first view considers that the magnitude of the stress required to nucleate twins is approximately the same throughout a sample despite the multiple sources and magnitude of local stress concentrations from which twins may nucleate. When the externally applied stress increases, there is an increasing number of stress concentrators (e.g., lattice defects) that reach some kind of critical activation stress so that twinning occurs. The alternative view considers that there are potential twin nucleation sites, each with different critical activation stress that reflects the strength of the lattice structure at that site. That means that twinning would not be simply determined by the macroscopically applied stress. However, if widespread twinning occurs over a small interval of increasing externally applied stress, it may be possible to define an apparent CRSS and estimate its value.
Other investigations aimed at empirically calibrating the CRSS for plastic twinning using inversion of calcite twins from experimentally deformed monophase natural samples (hence non-perfect polycrystals) under controlled stress conditions [29][30][31][32]. They have resulted in a relatively uniform value (10 ± 4 MPa) for the CRSS at a temperature below 150-200 • C, especially when grain size and strain hardening are taken into account ( Figure 4A). This value has been considered to be convincing by [31] and is widely used in calcite twin analyses [3,8,9,33]. Because the sources of stress concentrations such as grain-scale heterogeneities are very numerous in natural crystals (dislocations, fractures, indenters, pre-existing twins), this empirical CRSS may represent the critical stress required for the twin to propagate/grow when already nucleated at a stress concentrator and to become macroscopically visible. In other words, like in Griffith's theory of microfracturing, one can imagine nuclei as being there and that one just needs to activate them. The scatter in the measured CRSS value may be attributable to a range in the strength of the stress concentrations [34]. This is in line with the observations that a twin propagating across one grain may nucleate a twin in the adjacent grain when the local stress concentration associ-ated with its tip reaches the grain boundary, giving birth to bursts of twin development in aggregates.   [7,28,30]. (B) Stresses at the onset of elastic and permanent twinning in Carrara marble (grain size 150 µm). The curves show trends in the data and are not fits. In each experiment, the onset of permanent twinning was at an axial volumetric elastic strain of <0.1%, after [18]. (C) Evolution of the CRSS value as a function of grain size and strain. Modified and implemented after [13].
The above considerations would suggest that the empirical CRSS to be considered for the purpose of paleopiezometry would better be the stress associated with plastic twin propagation across the grain (controlled by external stresses and governed by a Schmid-type criterion [18,35]) rather than the stress associated with twin nucleation (likely controlled by local lattice structure and stress concentrations). As a result, the existence of a CRSS for twinning in calcite has more of an empirical than a theoretical basis (and a mechanistic interpretation as well). The notion that the CRSS is the stress at which a statistically significant number of nucleation sites become active, and that this is quite a sudden thing, provides a possible elegant way to reconcile the picture of the twin-forming process with the practical existence of a CRSS. Note that [13] question whether the CRSS value determined from the macroscopic stress-strain behaviour should be considered at the appearance of the first isolated plastic twin lamella (initiation) or at the "plateau" of the macroscopic stress-strain curve (nearly steady-state flow stress) which is marked by the densification and thickening of twin lamellae.
Following the earlier works by [28][29][30][31][32]36,37], experimental work over the last decade has been devoted to better constrain the stress conditions for twinning, including a better estimate of the value of the macroscopic (empirical) CRSS and its dependence on temperature, strain, strain-rate and grain size [13,18,23,38]. Figure 4 summarizes the information available to date about the possible value of an (apparent) CRSS for plastic twinning in calcite. Covey-Crump et al. [18] were the first authors to provide an estimate for the stress needed to initiate elastic twinning in calcite ( Figure 4B). The authors observed abrupt changes in the variation of strain with applied stress at some small stress (~10 MPa in Carrara at 200 • C) and at a greater stress (~65 MPa in Carrara at 200 • C, decreasing to 41 MPa at 500 • C). As the deformation was recoverable if the samples were unloaded when taken to a maximum stress less than 65 MPa, they thought in terms of elastic twinning commencing at~10 MPa and permanent twinning at~65 MPa ( Figure 4B). According to the authors, the reported value of up to~10 MPa is not a CRSS because the onset of elastic twinning within a grain is more strongly influenced by the number of its e-twin systems which are possibly activated (i.e., that have positive Schmid factors) than by the actual values of those Schmid factors. Furthermore, the value cannot be converted into a CRSS using the Schmid factor because it is not just a single twin system that is involved in any particular grain. In contrast, the onset of permanent twinning is governed by the value of the Schmid factor on the most favourably oriented e-twin system within the twinning grains, so assuming a most favourably oriented twin system with a Schmid factor of 0.5, a CRSS of~30 MPa can be inferred (at~200 • C for a mean grain size of 150 µm).
The contrasting differential stress values suggest that a permanent twin requires a different (and higher) CRSS than an elastic twin, the underlying reason being that elastic and plastic twinning involves a different accommodation of the deformation. Yet, it is still unclear if there is a transition from a thin elastic twin to a permanent twin lamella, i.e., if an elastic twin may be considered as some precursory development of a permanent twin. The uniaxial experiments performed by Parlangeau et al. [13] at room temperature on monocrystals with different grain sizes but overall much larger (~6 mm) than the mean grain size of natural aggregates yielded a very low CRSS, in the range 0.5-1.25 MPa ( Figure 4C).

Dependence of CRSS on Temperature and Strain Rate
The onset of elastic twinning should be significantly temperature sensitive, and the onset of permanent twinning may also well be if dislocation processes are important (Covey-Crump, personal communication, 2021). The observation that an apparent/empirical CRSS for the onset of mechanical twinning in calcite can be possibly defined and that it exhibits only a small decrease with increasing temperature ( Figure 4A) and decreasing strain-rate, matches the findings from deformation experiments on engineering materials [11].
Experiments show negligible changes in CRSS with strain rate over the range 10 −1 to 10 −7 s −1 [37]. The apparent CRSS for the onset of mechanical twinning in calcite, therefore, exhibits a limited dependence on temperature and strain rate [9]. Especially, this near independence of twinning upon strain rate is a definite prerequisite for the use of calcite twinning as a paleopiezometer.

Dependence of CRSS on Strain
The CRSS for twinning decreases from~12 ± 4 MPa at 20 • C to~6 ± 4 MPa at 400 • C [30,36,37] depending on the strain ( Figure 4A). The stress necessary for twin pro-duction in polycrystals, therefore, depends on total strain [9,28,30]. This is because calcite hardens once twinned: TEM observations show increased dislocation density near the twin boundaries. Repulsive interactions occur between twins and the dislocation structure. Twin lamellae hinder each other in their development. Such negative interactions result in strain hardening, during which both dislocation slip and twin growth are made more difficult (e.g., [23]).

Dependence of CRSS on Grain Size
Mechanical twinning is sensitive to grain size (e.g., [23,26,39,40]), but again this effect is poorly characterized [13]. Twinning in calcite is increasingly difficult as grain size decreases. Rowe and Rutter [26] have shown that while twin density (number of twins per mm) is poorly grain-size dependent, the twinning incidence (% of grains in a given size range containing optically visible twins) and the twin volume fraction (% of the volume of the twin fraction of the grain) are found to be largely dependent on grain size. Twinning is easier in large grains, which is related to the fact that grain boundaries are obstacles to the propagation and thickening of twins. As a result, with decreasing grain size greater stresses are required to ensure that the energy cost of creating new twin boundaries does not outweigh the mechanical work done by the twin [41]. Alternatively, finer-grained polycrystals require a greater twin number density to achieve a given strain, and this requires greater stress [42]. Given that the stress required for permanent twinning is linked to dislocation processes, a Hall-Petch-type grain size dependence can therefore be expected [39]: (σ 1 − σ 3 ) = σ 0 + kd −1/2 , where σ 0 and k are material parameters and d is grain size ( Figure 4C).
To evaluate the sensitivity of the stress required for twinning to grain size, Covey-Crump et al. [18] performed uniaxial compression experiments on a coarse-grained marble (Carrara marble, mean grain size of 150 µm) and a fine-grained limestone (Solnhofen limestone, mean grain size of 5 µm) mainly at 200 • C, i.e., at a temperature which does not exceed too much the range of temperatures encountered in the upper few kilometres of the shallow crust (e.g., sedimentary basins). The total range of investigated temperatures is 20 • -600 • C. Permanent twinning in the Carrara marble samples initiated under differential stress of 65 MPa at 200 • C, while twinning was not observed in the Solnhofen limestone samples deformed at 200 • C at applied stresses <180 MPa ( Figure 4B). This supports the dependence of twinning occurrence and CRSS on grain size.
Moreover, if twinning is sensitive to grain size and if twin propagation across one grain may nucleate a twin in the adjacent grain when the local stress concentration associated with its tip reaches the grain boundary, then twinning is probably also sensitive to grain size distribution. That would be because the number of grains in which twin nucleation can be easily stimulated will be influenced by the proportion of relatively small grains [43].
However, it is important to remember that the determination of the size of the analysed grains is usually made on a 2D surface whereas crystals are 3D objects. Unfortunately, there are only a few non-destructive techniques that can be used to routinely determine the grain size in 3D. Determining the mean size of grains in 2D, using for instance (length max + length min)/2, has a major drawback. The 2D appearance of a grain can be misleading since it could correspond to a small section of a larger grain. This issue may give rise to bias in the evaluation of the mean and range of grain sizes in natural samples, with a direct impact on the value of the CRSS to be adopted in stress calculation. Figure 4C compiles the available information about the possible evolution of the CRSS as a function of grain size, assuming the evolution of the CRSS against grain size follows the Hall-Petch rule and that CRSS also depends on strain/ internal grain deformation since calcite hardens once twinned [34,36,43].
Other possible issues, which remain to be further investigated, are related to the likely influence of the crystal chemistry as stresses necessary for twinning might also be affected by the Ca-Mg substitutions within calcite, the presence of second phases, nano-scale porosity, changes in pore-fluid chemistry, and solid solutes. It seems that at first glance, Ca-Mg substitutions within calcite do not influence twinning [44] but this needs to be confirmed by dedicated investigations in the future. Note however that the occurrence of brighter twins under cathodoluminescence in crystals displaying otherwise homogeneous luminescence ( Figure 1E) is unlikely to be related to pre-existing chemical heterogeneities or diagenesis. We propose that this observation reflects the local diffusion of activator elements that would concentrate during mechanical twinning. Albeit this hypothesis requires further chemical analyses to be validated, a possible candidate is Pb, as this element is often found in calcite and is a strong activator of luminescence [45]. Similar phenomena of Pb diffusion/concentration under local stress and strain have been described in other minerals [46,47]. Table 1 provides a qualitative summary of the estimated respective influence of differential stress, fluid and confining pressure, temperature, strain, strain rate, grain size and crystal chemistry on calcite twinning (twin morphology and development).

The Need for Accommodation of the Geometric Incompatibility at Twinned Grain Boundaries
Twinning leads not only to very inhomogeneous deformation inside the grain but also to geometrical complications at grain boundaries. Within a polycrystal, the ease for a grain to twin is therefore strongly dependent on the constraints exerted by the neighbouring grains. These constraints may include the ability of the surrounding grains to deform either elastically and/or plastically according to their crystallographic orientation, which will greatly condition the amount of shear that can be transmitted across the grain boundaries.
In response to 'simple shear' by twinning, and in absence of any recovery mechanism, twinned grains would expectedly exhibit angular steps at grain boundaries and impinge on neighbouring grains [4]. Such prominent step (or symmetrical corresponding trough) at the surface of a twinned grain is related to the twin lamellae reaching the surface as revealed by BEM ( Figure 5A,B). AFM (Atomic Force Microscopy) can be further used to measure the topography of the surface of the twinned grain, with the changes in height (slopes) of the surface being related to the occurrence of twin lamellae ( Figure 5C). The difference in height across the twin lamellae is of the order of few nanometres. In more detail, the peak force error diagram of Figure 5D even further shows the actual occurrence of several twin lamellae with variable thickness.
Some kind of accommodation must be achieved to prevent the formation of voids during deformation. In naturally deformed rocks, dissolution-recrystallisation along grain boundaries, governed by stress differences, help remove geometric incompatibilities and thereby prevent intense elastic deformation and the formation of open spaces and fractures. Twin boundary migration may tend to create rectified grain boundaries. Grain-boundary migration is a thermally activated process primarily driven by differences in stored energy between neighbouring grains and which could remove most of the geometrical incompatibilities [4]. At low temperatures, where twin thickening (by e.g., dislocation climb) may be usually limited, increasing shear strain will be accommodated by nucleation of new twins instead of broadening of existing twins, resulting in a high density of thin twins. These thin twins are presumably associated with steps at grain boundaries which are significantly smaller than those related to the thicker twins which are expected to form at higher temperatures [31]. At low temperature, pressure-solution is therefore the most likely operating mechanism that would remove such steps and ensure geometrical compatibilities at grain boundaries.
Parlangeau et al. [13] reported that their experimental results on the formation of thick twins in unconfined grains deformed at room temperature cannot be easily extrapolated to natural aggregates. This is because confining the samples with a fluid or leaving an unconfined free surface as during uniaxial compression could affect twin thickening with respect to confinement by neighbouring grains. At low temperatures, shear strain compatibility at grain boundaries would expectedly favour the occurrence of numerous distributed thin twins. Conversely, the presence of fluid confined or free surfaces could favour the development of thick twins.

Calcite Twins: A Low-Temperature Geothermometer?
The thickness of twin lamellae is often thought to be mainly controlled by temperature [4,48]. Four types of calcite twins can be recognized: type I twins are thin (<1 µm), straight and dominate at temperatures <170 to 200 • C; type II twins are thicker (>1-5 µm) and slightly lens-shaped, and prevail at temperatures between 150 and 300 • C; type III twins are thick, curved and tapered and form at temperatures >200 • C; type IV twins are thick, patchy, irregular, accompanied by grain boundary migration and occur at temperatures >250 • C. A given amount of strain is expected to be accommodated by a number of thin twins at a temperature lower than 170-200 • C and by fewer but thicker twins at a temperature above this range. This has led to the proposal of a calcite twin-based low-temperature paleothermometer [5] (Figure 6), which received great success in allowing to roughly assess the temperature at which twin strain occurred (but obviously not allowing to assess the temperature of calcite precipitation since the timing of deformation within the crystal lifetime is unknown).
ture [4,48]. Four types of calcite twins can be recognized: type I twins are thin (<1 µm), straight and dominate at temperatures <170 to 200 °C; type II twins are thicker (>1-5 µm) and slightly lens-shaped, and prevail at temperatures between 150 and 300 °C; type III twins are thick, curved and tapered and form at temperatures >200 °C; type IV twins are thick, patchy, irregular, accompanied by grain boundary migration and occur at temperatures >250 °C. A given amount of strain is expected to be accommodated by a number of thin twins at a temperature lower than 170-200 °C and by fewer but thicker twins at a temperature above this range. This has led to the proposal of a calcite twin-based lowtemperature paleothermometer [5] (Figure 6), which received great success in allowing to roughly assess the temperature at which twin strain occurred (but obviously not allowing to assess the temperature of calcite precipitation since the timing of deformation within the crystal lifetime is unknown). Figure 6. Twin strain vs. twin width plot illustrating the temperature dependence of twin geometry that may be used as a low-temperature paleothermometer. Modified after [9,31].
The tendency for calcite grains to develop thick twins at a temperature over 200-250 °C can be explained by considering the physical mechanisms of inelastic deformation. Twin thickening within a grain should be assisted by recovery mechanisms, such as dislocation recovery, which would reduce repulsive interactions between twin dislocations and matrix dislocations [23]. Dislocation recovery mechanisms, grain boundary migration, grain boundary sliding or even pressure-solution relieve strain incompatibilities at grain boundaries (see Section 2.3; [4]; Figure 5). These mechanisms are favoured by high temperature and low strain rates. As temperature increases, the hardening coefficient is reduced and twin thickening prevails [23].
Rowe and Rutter [26] suggested that the average twin width increases linearly with grain size. However, although temperature might be the major influence determining twin width, average width might be also influenced by changes in strain (and duration of application of stress as well), strain rate, or grain size.
Janssen et al. [49] reported that twin morphology and twin width do not correspond to deformation temperatures at least for temperatures below 250 °C and that deformation Figure 6. Twin strain vs. twin width plot illustrating the temperature dependence of twin geometry that may be used as a low-temperature paleothermometer. Modified after [9,31].
The tendency for calcite grains to develop thick twins at a temperature over 200-250 • C can be explained by considering the physical mechanisms of inelastic deformation. Twin thickening within a grain should be assisted by recovery mechanisms, such as dislocation recovery, which would reduce repulsive interactions between twin dislocations and matrix dislocations [23]. Dislocation recovery mechanisms, grain boundary migration, grain boundary sliding or even pressure-solution relieve strain incompatibilities at grain boundaries (see Section 2.3; [4]; Figure 5). These mechanisms are favoured by high temperature and low strain rates. As temperature increases, the hardening coefficient is reduced and twin thickening prevails [23].
Rowe and Rutter [26] suggested that the average twin width increases linearly with grain size. However, although temperature might be the major influence determining twin width, average width might be also influenced by changes in strain (and duration of application of stress as well), strain rate, or grain size.
Janssen et al. [49] reported that twin morphology and twin width do not correspond to deformation temperatures at least for temperatures below 250 • C and that deformation temperatures that are derived from twin width and twin morphology are not in good agreement with temperatures estimated by other independent methods (e.g., vitrinite reflection, conodont colour alteration index). It comes that although the width of calcite twin lamellae does change as a function of deformation temperature, this change is seemingly not systematic.
Curzi et al. [50] reported field data from the Mt. Tancia Thrust in the Central Apennines. They showed that during compressional deformation, the thrust tip progressively propagated upward and that the shear zone was associated with a pervasive fabric that developed at temperatures of 50-70 • C as indicated by clumped isotopes from associated veins and in agreement with the maximum burial depth of the investigated sedimentary units. This study reveals a discrepancy between geothermometers since the ambient temperature was significantly lower than the 150-200 • C previously inferred from the type II morphology of calcite twin lamellae.
The experiments by [23] further questioned the statement that the thick twinning occurs at a temperature higher than 170 • C by documenting few µm-thick twin lamellae formed during their experiments carried out at room temperature. These authors also show that the thickness of twin planes can vary within the same family of twin planes. The mean width of calcite twins increases with both temperature and strain, and thus, measurement of twin width provides only a rough estimate of peak temperature, unless additional constraints on deformation are known [23]. This finding has been later confirmed by the experiments on single grains by [13] showing that the duration of stress application, and hence the resulting strain, has also a significant impact on twin lamella thickness even at low temperature. Note, however, that in the latter study a large amount of thickening observed during the experiments is likely due to the lack of confining pressure, hence to the free surfaces of the single crystals.
Relating thoughtfully and reliably the occurrence of thick twins to increasing temperature would require detailed knowledge of the elementary processes involved in twin nucleation, twin growth, and the interactions of matrix dislocations with twin boundaries and twinning dislocations. Especially, the kinetics of these processes is required to predict the strength at very slow strain rates. If twin nucleation and twin growth have distinctly different kinetics, then the temperature at which thick twins prevail will also depend on strain rate [23]. This needs to be taken into account for the interpretation of natural twin data.
Schuster et al. [15] carried out experimental deformation of polycrystalline calcite to high strain at room temperature and confining pressures of 1-4 GPa using high-pressure torsion. The authors report that the experimentally produced mechanical e-twins are significantly thicker than 1 µm and have a moderate lenticular shape. This morphology rather conforms to morphologies typical for substantially higher deformation temperatures, e.g., type II twins [4]. The formation of thick, lens-shaped twins requires that the dislocations that propagate the twin can move out of their glide plane by climb or cross-slip [12]. Although this process is usually considered to be thermally activated, the results by [15] show that high lattice strain produced during the deformation under a confining pressure in excess of 1 GPa can lead to substantial twin thickening even at room temperature. Thus, although such confining pressure conditions largely exceed those under which deformation twinning considered in this paper may occur, this effect is possibly not as negligible as previously thought [51] and should not be neglected when calcite twin morphologies are used as geothermometers.
As a result, these contradictory observations plead in favour of more caution in using calcite twin morphology only (e.g., thickness) as a straightforward and precise lowtemperature deformation paleothermometer. Whereas thin twins appear to be diagnostic of low temperature (<170 • C) deformation, thick twins do not always reflect higher temperature (>170-200 • C) deformation, so twin morphology should not be blindly used as a valid paleothermometer, even when strain constraints are available.

Calcite Twin Data Measurements for Paleostress/Paleostrain Analyses
In the diagenetic domain of sedimentary basins, calcite grains in veins and host rocks may be large (100-500 µm) and may contain very thin twins (thickness < 1 µm) owing to weak twinning strain (<2-4%) achieved at temperature <150-200 • C ( Figure 1). Twin data are commonly collected optically from two or three perpendicular thin sections using a U-stage. 25-30 calcite grains are generally measured in each thin section [52]. Depending on the purpose of the analysis, the orientations of the C-axis and the three potential e-twin planes are measured for each crystal and the twinned/untwinned status of each twin plane is evaluated. The calcite twinning strain gauge technique [53] and the oldest techniques for stress analysis [54] only rely on optically visible twins, while more sophisticated inversion techniques also rely on the orientations of untwinned planes. The untwinned potential twin planes, along which twinning was not activated, hence which are associated with no optically visible twin lamella, can obviously not be directly measured, but their orientation can be deduced from the known symmetries and angular relationships between crystallographic planes and C-axes in calcite crystals. The measurement of grains without any or with a single twinned plane is therefore feasible if cleavage planes are present, but impossible otherwise.

Universal Stage
The classical technique for measuring calcite twins is the use of a U-stage. The accuracy of U-stage measurements of microstructures of a grain is usually about 1-4 • [55], and the accuracy of measuring directly a pair of calcite host C-axis and the related twin plane is about 9 • [56]. The direct measurement of the C-axis under the U-stage is, however, timeconsuming and its determination may be inaccurate, e.g., in the case of pseudo-undulating extinction. In some instances, measurements can be assisted by a software that significantly speeds up the measurements and makes them more accurate. The software calculates the orientation of the C-axis from the orientations of the optically visible crystallographic planes (twin or cleavage planes) using the angular relationships between the C-axis and the crystallographic planes in the crystal [57,58]. While bringing true progress in data collection, the software [57] does not allow one to collect data from a crystal having only one visible plane (twin or cleavage). This software (and the TwinCalc package as well) check the angles between C-axes and poles to e-twin planes during input and notifies the user if the difference between the theoretical and the measured angles exceeds 3 • . To the latter value must be added the standard usual optical error which is estimated about 1-2 • .
The tilt of a thin section mounted on the U-stage is limited by the bumping of the U-stage into the objective. This maximum tilt is generally 50-60 • (e.g., [59]). Moreover, the thin section should even not be tilted at this maximal angle because the error on the orientation measurement caused by refraction increases with increasing tilt. As a result, the U-stage only allows for the safe measurement of a relatively small area of the thin section and moderate dipping and plunging microstructures (~<50 • ) relative to the thin section. This limitation, therefore, affects the number of twin planes for which the orientation can be accurately measured and the twinned or untwinned character can be unambiguously checked. In addition to being long and tedious, data collection by optical means, therefore, suffers from the need for ensuring that potential twin planes are twinned or not, which may be difficult to achieve for some of them. This has a direct implication for the accuracy of the inversion of calcite twin data for stress.
The percentage of potential twin planes oblique at a too low angle to the thin section plane to be reliably measured, hence possibly misclassified as untwinned planes (i.e., biased) has been estimated by [60]. For a virtual horizontal thin section with 100 randomly oriented grains (i.e., 300 e-twin planes), and considering that a twin plane lying at an angle of 30 • or less to the thin section cannot be measured using a classical U-Stage, the percentage of potentially optically biased twin planes is estimated to be lower than 11%, with an average of 6%. Note however that optical bias may also result from the number, density, width and spacing of twin lamellae as well as the degree of extent of lamellae across the grains, which may result in an even greater percentage of possibly biased twin data (up to 25%, [61]).
The use of U-stage measurements in petrofabric analysis, including measurements of calcite twins, has been recently revisited [62].

Automated Fabric Analyser
Köpping et al. [63] reported the use of a fabric analyser to speed up and improve the collection of calcite twin data to be used in further paleostress and paleostrain analyses. A fabric analyser is an automated polarizing optical microscope that allows for the determination of the orientation of the C-axis as well as the twinned part of a calcite grain at each pixel in the field of view. The output of the instrument is a set of images and a pixel data file for every single measurement [64]. The fabric analyser measures a complete thin section with~5 µm/pixel resolution in less than one hour and without any specific sample preparation. Measurements of C-axes and twin lamellae orientations in calcite grains as collected with the fabric analyser are used as input in a software package (named PACT: Paleostress Analysis with Calcite Twins) to determine both paleostress (using Turner's technique) and paleostrain (using Groshong's technique) from twin data.
The authors conclude that their so-called 'automated paleostress analysis' is a significant step forward compared to the use of a U-stage. This conclusion seems to be oversold. Clearly, this device has two advantages: (1) there is no limitation in measuring C-axes orientations even for crystals with plunges ≥50 • to the analysed thin section, and (2) the use of the fabric analyser saves time for the collection of visible twins and related C-axes in calcite grains since [63] report that 50 host/twin pairs can be measured in about 30 min. Depending on the grain size and the twin thickness, such minimum measurements may well be sufficient for a rough paleostress analysis consisting simply of the use of the Turner's technique, with maxima of the statistical orientation of P and T-axes estimated for each C-axis-twin lamella pair and interpreted as possible locations of σ 1 and σ 3 [56] or for Groshong's analysis. However, the gain in time is limited, and even null if one aims at also measuring the orientations of untwinned planes (along which by definition no twin lamella can be observed) as it is the case of more sophisticated calcite twin inversion techniques for stress [60,[65][66][67]. Finally, as stated by the authors, one constraint of the fabric analyser is the maximum resolution of~5 µm per pixel which bounds the measurable twins to a thickness of a minimum of 10 µm or 2 pixels, correspondingly. This means that the device can seemingly not allow for measurements of thin twins as commonly encompassed in limestones deformed at a temperature lower than 170 • C (e.g., reservoirs).

EBSD
The EBSD (Electron Back-Scatter Diffraction) is commonly used to obtain quantitative information about microstructures and crystallographic orientations of mineral phases, hence extensively investigate the fabrics of metals and rocks displaying very small grains. EBSD reveals the grain size, grain boundaries, crystallographic orientation and texture of the mineral phase. EBSD protocol requires a specific sample surface polishing step but makes it possible to work on a 3D sample and not necessarily on a thin section. Each face of the sample can be analysed separately. The orientations of the twin planes can be measured regardless of the dip with respect to the analysis surface and with an accuracy of 1 • .
EBSD can theoretically help solve the U-stage measurement issues, hence improve the efficiency of data collection of twinned and untwinned planes. However, because calcite grains in veins and sedimentary host rocks may be large (100-500 µm) and may contain very thin twins (thickness < 1 µm), twin data acquisition using EBSD faces the problem of the number and spacing of required spots against the low thickness of twins and the number and large size of the grains to be measured ( Figure 7A).
To date, only a few authors have reported the use of EBSD to collect calcite twin data for stress inversion [68][69][70][71] and measurements were taken from a limited number of grains. Rez et al. [68] reported that a complete EBSD analysis of a 1 cm × 2 cm sample would take weeks to obtain usable results for inversion. They, therefore, propose to use an analysis grid composed of perpendicular lines extending over the entire surface of the sample with a specific spacing and measurement step between successive spots. Evenson [69] used EBSD to determine the crystallographic orientations of grains and also to determine the thickness of the twin lamellae. Because the twins were very thin (<1 µm width), a step varying from 0.1 to 0.5 µm was used.   (1), the principle of twin data collection along perpendicular lines in grains with a step of 0.5 µm between successive measurement spots (2) and the principle of the identification of the slip system (here, e-twinning) (3,4). (B) Compared pros and cons of twin data collection using U-stage and EBSD.
Obtaining the information regarding crystal orientation on a single crystal of about 100-300 µm in diameter would therefore require on average several days using a 0.5 µm step. If only portions of crystals are analysed or depending on the measurement step, there is a risk for a family of twins to be missed. The main problem with EBSD is therefore the measurement time.
Parlangeau et al. [72,73] investigated a new EBSD-based procedure to collect twin data, which is a compromise between improving twin data measurement quality compared to U-stage and minimizing the time of data acquisition. For a tectonic vein, it consists of collecting twin data along perpendicular lines in grains-which requires the preliminary identification of the grain boundaries-with a step of 0.5 µm between successive measurement spots ( Figure 7A). The efficiency of this new procedure was tested against U-stage measurements on the same grains. Three cases were encountered, with EBSD providing the same, better or worse twin detection/measurements than U-Stage. The latter case is related to the missing of some very thin or heterogeneously distributed twin lamellae, due either to the too large EBSD measurement step or to the strategy of data collection along lines. Although this case is likely to be rarely met, it is impossible to exactly determine how frequently it may occur in a given thin section, which casts some doubts on the reliability of this procedure.
As a result, after comparing EBSD-based and U-Stage-based measurement protocols, one cannot conclude about the real advantage of EBSD ( Figure 7B). This is because EBSD analysis is entirely surface dependent and only explores the sample to a thickness of a few nanometres. That is, if a plane does not cross the entire crystal and ends up on the surface analysed by EBSD, it will not be detected as an activated twin plane. The advantage of the U-stage over EBSD is the ability to see and therefore measure activated crystallographic plane by transparency. Another issue of not scanning the entire grain by EBSD consists in missing a twin which would not extend across the grain. Making measurements along lines and not scanning the entire surface of the grain can therefore lead to misevaluate the twinned status of twin planes. This adds to the debatable choice of 0.5 µm measurement step which may be too large for the reliable detection of thin twin lamellae (thickness < 0.5 µm). If time acquisition is not an issue, however, EBSD reveals to be the best technique to collect twin data provided that the whole surface of the grains is automatically scanned with a measurement step not larger than 0.1 µm.
Even if calcite twin data collected using EBSD were not used for stress inversion purposes, McNamara et al. [74] recently combined EBSD and chemical mapping of calcite veins to characterize calcite vein fill morphologies and to gain insight into the mechanisms of calcite fracture sealing and reservoir's history of stress, strain, and deformation in the Kawerau geothermal field in the basement of the Taupo Volcanic Zone, New Zealand.

Calcite Twins: Indicators of Stress, Strain or Both?
Stress and strain are different physical entities. Whether calcite twins can be treated as strain or stress indicators is a long-lasting debate, recently revived by [61]. Paleostress analysis using calcite twins received criticism in that twinning inherently represents (finite) strain (e.g., [75,76]), i.e., geometrically analogue to simple shear of grains. Strain analysis obviously relies on optically measurable twins. It is the same for the old techniques of stress analysis (e.g., [54,77]), while more recent techniques take also into account untwinned planes (along which by definition no twinning deformation has occurred, hence is not visible). Methods of strain/stress analysis based on calcite twins share the fundamental assumption that the measured twins formed in a homogeneous stress field and were not passively rotated after formation. Burkhard [4] emphasized that these methods are best applied to very small strains that can be approximated by coaxial conditions so that the orientation of small twinning strain can be reliably correlated with paleostress orientation [76].
As emphasized by [61], both strain and stress approaches are considered valid on an empirical basis only, so none of them has a more robust theoretical support of its ability to yield correct results. The orientations of twinning strain tensor axes are empirically considered to reflect macroscopic tectonic deformations [78]. The orientations of twinning stress tensor axes are empirically considered as reflecting regional stresses based on their compatibility with other indicators (faults and fractures) for which the relevance to the regional case has been established for long (e.g., [52,79,80]). Strain and stress fields in polycrystalline calcite are heterogeneous [4,18,23,81] but Spiers [81] has claimed that deformation is distributed inhomogeneously between grains while the stress is much more homogeneous at the scale of the aggregate. So none of the approaches surpasses the other and is better theoretically supported than the other.
One can wonder whether the stress and strain determined from e-twins in a sample of polycrystalline calcite do have similar orientations of principal axes and similar shape ratios, especially in the case of polyphase tectonics with different stress conditions, even if the low strain hypothesis is fulfilled. Amrouch et al. [82] were the first authors to carry out separate stress and strain analyses from the same set of twin data from veins and host rocks in the Sheep Mountain anticline using CSIT (see Section 6) and CSGT (see Section 8). In terms of orientations of shortening ε min /compression σ 1 axes, the results show very good consistency, for both layer-parallel and late fold tightening preceding and following immediately fold growth, respectively, therefore supporting that internal strain of folded strata remained mainly coaxial [82] although the stress/strain ratios display some differences. Wakamori and Yamaji [83] more recently tackled the question of strainstress coaxiality by combining the [84] strain analysis method with a novel method of paleostress analysis of mechanical twins, which clusters the directional data of e-twins using a statistical mixture model and determines stresses for each group of data. They applied this integrated stress-strain analysis to calcite veins from a Miocene basin of central Japan and like [82], they obtained close principal stress and strain axes (<20 • ) as well as similar values of the stress/strain ratios. It is however not straightforward to generalize this result, which should be theoretically checked case by case.
Even researchers who have for years reported calcite twin strain analysis now consider that although the result is actually a strain tensor, a similar orientation of the stress tensor appears likely [85] if deformation remains small and/or irrotational (coaxial). However, one can wonder whether such inference is sound in the vicinity of, or inside, fault zones where the analysis has been recently carried out [85] (see Section 8).

Use of Calcite Twins as Markers of Stress Regime: Inversion of Calcite Twin Data for Principal Stress Orientations and Stress Ratio
Among the available techniques that analyse calcite twins in term of stress, some are designed to determine paleostress orientations only [54,77,86], paleostress magnitudes only [26,87] or both [58,60,[65][66][67]88,89]. In this section, we consider only the techniques using calcite twins to determine principal stress orientations (and stress ratios). The techniques using calcite twins as stress gauge, i.e., to determine differential stress, will be considered specifically in Section 9.

Brief Review of Recent Methods and Improvements of Calcite Twin Inversion for Stress
While Burkhard [4] simply listed the theoretical pros and cons of the twin analysis methods available at that time [54,65,77,88,89], Shelley [90] was the first author to provide a true comparison of three available methods. Later, others [3,9,66,[91][92][93] emphasized the respective limitations of previous methods. It is only recently that a real effort has been dedicated to proposing new methods/inversion schemes to go beyond older methods [58,60,66,67]. To the best of the authors' knowledge, however, the technique by [66] has to date been applied neither to experimentally or naturally deformed samples nor in extensive regional studies. After recalling the principle of the widely used CSIT, the new techniques set out during the last decade are briefly described hereinafter.

Etchecopar Technique: CSIT
The inversion process of CSIT is very similar to that used for fault-slip data [65] since twin gliding along the twinning direction within the twin plane is geometrically comparable with slip along a slickenside lineation within a fault plane. The method relies on the existence of a critical resolved shear stress (CRSS) for twinning (see Section 2.2) and on the assumption that a potential twin plane is twinned (respectively, untwinned) if the resolved shear stress applied on it is greater (respectively, lower) than the CRSS. In contrast to most earlier methods, this inversion process takes into account both the twinned and untwinned planes.
The principle of the inversion of calcite twin dataset is to find a stress tensor (or several stress tensors) which verifies these two inequalities for the largest number of twinned planes and the whole set of untwinned planes. The solution has the form of a reduced stress tensor with 4 parameters: the orientations of principal stress axes (σ 1 , σ 2 , σ 3 ) and the stress ratio [=(σ 2 − σ 3 )/(σ 1 − σ 3 )]. The differential stress (σ 1 − σ 3 ) is normalized to 1, so the normalized resolved shear stress applied on each twin plane varies within The inverse problem consists of finding the stress tensor that best explains the orientation distribution of measured twinned and untwinned planes. The first step consists of an arbitrary choice of a percentage of twinned planes to be explained. Then there is a random trial of a number of reduced stress tensors. The resolved shear stresses are calculated on the twinned and untwinned planes, which are ranked as a function of the decreasing resolved shear stress. The solution tensor should theoretically meet the requirement that all the twinned and untwinned planes should be consistent with it. Thus, all twinned planes should sustain a resolved shear stress (τ s ) larger than that exerted on all the untwinned planes. The sorting allows one to determine rapidly whether some untwinned planes are incompatible with the tensor (i.e., the resulting resolved shear stress is greater than that for some compatible twinned planes). The fit of the stress tensor to the dataset is evaluated through the use of a penalization function, f, ideally equal to 0, which is defined as: where τ s min is the smallest resolved shear stress applied on the twinned planes compatible with the tensor, and τ s j is the resolved shear stresses applied on the j untwinned planes such that τ s j > τ s min . The penalization function increases if incompatible untwinned planes are incorporated in the solution. The optimal tensor is obtained when the maximum number of twinned planes and the minimum number of incompatible untwinned planes are incorporated in the solution.
This process, therefore, yields the orientation of the 3 principal stress axes, the stress ratio and non-dimensional differential stress, (σ 1 − σ 3 )/τ a = 1/τ s min . Under the assumption of a known and constant CRSS, τ a , the actual differential stress is given by (σ 1 − σ 3 ) = τ a /τ s min , so the technique provides the 5 parameters of the deviatoric stress tensor, i.e., principal stress orientations and differential stress magnitudes. The application of the CSIT technique to experimentally deformed natural samples [29,30] has returned maximum deviations of principal stress orientations compared to the experimentally applied ones of 5-7 • for monophase and 7-11 • for polyphase cases.
In contrast to its derivatives (e.g., CSIT-2, see below), the sorting of successive stress tensors from a heterogeneous (polyphase) dataset is done in sequence. If more than~30% of twinned planes in a sample are not explained by a unique stress tensor, the inversion process is repeated with the uncorrelated twinned planes and the whole set of untwinned planes to try to identify a second stress tensor. Because this procedure may lead to a loss of information since the removed twins consistent with the first stress tensor may well be consistent with the second tensor, Rocher et al. [94] modified the method by testing the twinned planes consistent with the first tensor against the second tensor: if some twinned planes are also consistent with the second tensor then they are incorporated into the second subset and the calculation of the second tensor is refined with the new dataset.

Yamaji Technique
Yamaji [61] has renewed the approach of inversion of calcite twin data for stress. The preliminary exploration of the extent to which calcite twinning may constrain stress [61] demonstrates that twinned planes better constrain the stress tensor than untwinned planes and, unexpectedly, that differential stress estimates are poorly resolved for differential stresses >50-100 MPa, a result later confirmed by [60,67].
A significant contribution of [61] has been to provide an elegant way for conceptualization and visualization of the stress solution and related consistent twin datasets. Based on the five-dimensional stress space that fulfils the principle of coordinate invariance, [61] shows that the inversion of twinned and untwinned planes is comparable with fitting a spherical cap to data points on a unit sphere S in the space. The principal stress orientations and the stress ratio (i.e., the reduced stress tensor) are indicated by the centre of the cap, whereas differential stress is denoted by the size of the cap, defining the deviatoric stress tensor responsible for twinning. Superimposed stress tensors have each their corresponding spherical cap, with a possible non-zero intersection if some twinned planes are consistent with at least two distinct stress tensors. The extension, hence the radius of the caps is constrained by the untwinned planes ( Figure 8).
Based on this geometrical interpretation Yamaji [66] proposed a new inversion technique based on the generalized Hough transform to determine the stresses from heterogeneous calcite twin data. The method attempts to fit a spherical cap to the distribution of compatible twin planes and simultaneously attempts to fit the complementary region of the cap on S to the distribution of compatible untwinned planes for a given stress. Like CSIT the 5 parameters of the deviatoric stress tensor are determined. The method is demonstrated with synthetic data sets to be robust to sampling bias (i.e., twinned planes misclassified into untwinned planes), variability of the CRSS and heterogeneity of data. According to [66], the treatment of deviatoric stress tensors in a five-dimensional space circumvents some possible issues of the earlier techniques [65,88] and makes the separation of superimposed stress tensors easier.
The present method separates stresses from heterogeneous data one by one by counting out twin data compatible with previously separated stresses. But since a twin datum can be explained by many stresses, the counting out reduces the information contents of the dataset so that the second solution becomes inaccurate. The method assumes that the spherical caps representing the stresses to be detected have no intersection, giving rise to the decreasing number of twin data and the diminishing accuracy during the detection of the second and following solutions. This is illustrated by the inaccuracy of the second solution in the most complex case tested (superimposed normal and strike-slip faulting stress regimes sharing a similar orientation of σ 3 and a stress ratio of 0.5 while the differential stress related to the first stress regime is twice than that of the second stress regime). This technique therefore efficiently separates superimposed stress tensors from heterogeneous datasets as long as the spherical caps corresponding to the stresses to be detected have no or a small intersection. The basic principle of CSIT-2 is similar to that of CSIT in that both techniques rely on the same algorithm of mitigation of the penalization function, but circumvents major limitations of CSIT. The major differences lie (1) in the systematic scanning of the space with a regular interval of 10° for the 3 Euler's angles (defining the orientation of the principal axes of the stress tensor) and a fixed stress ratio value of 0.5 instead of a random trial of a limited number of reduced stress tensors, and (2) the automatic detection of one or several tensors with similarly low penalization function 'in parallel', in contrast to the 'in sequence' identification of stress tensors with CSIT in which the first solution tensor retained for further optimization is one among the many possible tensors with the same low penalization [60]. The selection of tensors passes through an optimization process [106] to tweak the parameters of the stress tensor (orientation of principal stress axes and stress ratio). The ability of the new technique to detect, separate and determine stress tensors from monophase and polyphase twin datasets including measurements errors or various grain sizes has been demonstrated by testing it on both numerically generated twin datasets as well as naturally deformed polyphase samples [60]. The usual maximum methodological uncertainty for stress determination associated from a naturally deformed sample with a homogeneous grain size using CSIT-2 has been estimated ±10° (at very worst ±15°) for principal stress orientations and ±0.1 (=10%; at very worst ±0.4 (=40%)) for the stress ratio. For differential stresses it can be estimated ±37%; this rather high value is The basic principle of CSIT-2 is similar to that of CSIT in that both techniques rely on the same algorithm of mitigation of the penalization function, but circumvents major limitations of CSIT. The major differences lie (1) in the systematic scanning of the space with a regular interval of 10 • for the 3 Euler's angles (defining the orientation of the principal axes of the stress tensor) and a fixed stress ratio value of 0.5 instead of a random trial of a limited number of reduced stress tensors, and (2) the automatic detection of one or several tensors with similarly low penalization function 'in parallel', in contrast to the 'in sequence' identification of stress tensors with CSIT in which the first solution tensor retained for further optimization is one among the many possible tensors with the same low penalization [60]. The selection of tensors passes through an optimization process [106] to tweak the parameters of the stress tensor (orientation of principal stress axes and stress ratio). The ability of the new technique to detect, separate and determine stress tensors from monophase and polyphase twin datasets including measurements errors or various grain sizes has been demonstrated by testing it on both numerically generated twin datasets as well as naturally deformed polyphase samples [60]. The usual maximum methodological uncertainty for stress determination associated from a naturally deformed sample with a homogeneous grain size using CSIT-2 has been estimated ±10 • (at very worst ±15 • ) for principal stress orientations and ±0.1 (=10%; at very worst ±0.4 (=40%)) for the stress ratio.
For differential stresses it can be estimated ±37%; this rather high value is expected for high differential stresses. The usual uncertainty is expected to be about ±25-30% [60].

Shan et al. Technique
The first recent contribution by Shan et al. [86] was to introduce the concept of synthetic slip plane and to apply it in the dynamic analysis of calcite twin data. The synthetic slip plane is the linear combination of a pair of twinned and untwinned e-planes in a single calcite crystal. It is independent of either twinned e-plane or untwinned e-plane, and it has a similar shear stress sign to the twinned e-plane. As such, this auxiliary slip plane plays nearly the same role in graphic and numerical dynamic analyses as the twinned e-plane does. This auxiliary slip plane can be used together with the twinned e-plane either to further constrain the distribution of compression and tension axes in graphical dynamic analysis or to better estimate the reduced stress in numerical dynamic analysis. Meanwhile, a novel method is proposed to estimate the reduced stress tensor by maximizing the sum of the shear stress along the slip direction on each individual twinned e-plane or synthetic slip plane. The first application of this method to a series of numerically generated datasets has provided generally accurate and robust stress results, but one could wonder about the further use of this approach knowing the availability of more complete-although more complex-stress inversion techniques.
Shan et al. [67] inversion method for polyphase twin data is a compromise between CSIT and the method by [89] which considers minimizing the sum of the resolved shear stresses on incompatible e-planes on one hand, and the methods by [66,93] which aim at finding the solution of an optimization problem that considers not only incompatible twinned and untwinned e-planes but also compatible twinned and untwinned e-planes on the other hand. The method also calculates deviatoric stress tensors (reduced stress tensors and differential stress magnitudes) but differs from CSIT and CSIT-2 in that it separates related data by maximizing the sum of the resolved shear stresses along the compatible twinned e-planes under each reduced stress while minimizing the number of the incompatible (twinned and untwinned) e-planes shared by the unknown stresses. To facilitate the search of the reduced stress tensor in 5-dimensional sigma space, a constant CRSS value, 10 MPa, is used to rescale the tensor, but the real deviatoric stress tensor can be retrieved with the knowledge of the real CRSS. The authors adopted the evolutionary algorithm to solve the above optimization problem. Application of this method to a series of two-phase artificial twin datasets indicates the relatively high accuracy of stress estimate. One issue of this method, as stated by the authors, is that it takes an enormous amount of time to calculate. The method has been refined by also maximizing the number of compatible untwinned e-planes shared by the unknown stresses [107].

Rez Technique: TwinCalc
TwinCalc [58] is a freeware application for stress analysis based on calcite twinning. On the one hand, TwinCalc provides a number of tools for (1) twin data acquisition using an optical microscope equipped with a U-stage; (2) statistical analysis and plotting of twin data; and (3) stress analysis of twin data using most of the classical methods ( [54,108] for stress orientations and [26,87] for differential stress). On the other hand, CSIT has also been implemented in the software for the determination of stress orientations and differential stress magnitudes, including the most recent modifications of this method by [60,93], and some new modifications further enhancing the resolution of the method. A tool for stress phases similarity analysis based on a cluster analysis of 9D stress-tensor-vectors angles [109] is also available as a part of the application.
Rez and Melichar [93] proposed some improvement of CSIT by introducing a systematic search of a large number of reduced stress tensors instead of choosing a 'best-fit' reduced stress tensor from a number of randomly generated ones. Rez and Melichar [93] further claimed that the initial CSIT penalization function has several issues, especially rather blunt minima and producing a lot of "noise" and scattered false minima, with implications for the accuracy of tensor separation from polyphase datasets. For [93], the penalization function of CSIT is strongly dependent on the compatible twinned and the incompatible untwinned planes with the stress tensor solution. In addition, the space of solutions with a penalization function of 0 is too large. So, they proposed a new penalization function with sharper maxima, depending on the number of compatible twinned planes, the number of compatible untwinned planes and the number of incompatible untwinned planes.
The resolution of this method would theoretically be further improved by using τ b (defined as the highest resolved shear stress of all untwinned e-planes) instead of τ a (=τ s min , the smallest resolved shear stress along all compatible twinned e-planes, used in CSIT and CSIT-2 as the twinning threshold [58]). In a perfect sample, both values would be such as τ b ≤ τ a . However, stress analysis of natural data, which are rarely homogeneous, often yields some incompatible untwinned e-planes. Occurrence of an untwinned plane receiving a high resolved shear stress (i.e., well oriented for twinning) can be related either to grain-scale stress shadow or heterogeneity within the aggregates or to an error during data acquisition (biased twinned plane = misclassified e-twinned plane as untwinned e-plane). Occurrence of untwinned e-planes with resolved shear stress low but higher than τ a occurs as a result of the tendency to maximize the number of twins consistent with a given stress tensor. As emphasized by [58], in a polyphase dataset, some twinned e-planes may be directionally consistent with two stress tensors, but should not be assigned to both, possibly instead only to the stress tensor for which the twin receives the highest resolved shear stress. If one assigns as many e-lamellae systems as possible to one stress phase, there is a risk that to do so, twins with low resolved shear stress are incorporated, which lowers τ a hence artificially increases the differential stress, while some incompatible untwinned e-planes are incorporated in the solution while they should theoretically not be. Using τ b as the normalized CRSS instead of τ a would theoretically eliminate incompatible untwinned planes, and therefore increase the numerical quality of the tensor and yield smaller differential stress. As stated by [58], the real normalized threshold for twinning is probably somewhere between τ b and τ a , however numerically undetectable. However, τ b is obviously very sensitive to biased twinned e-planes, so a single misclassified untwinned e-plane receiving a high resolved shear stress can artificially shift τ b toward a high value. As a result, it is of prime importance to remove from the dataset all erroneous data, but since some (all?) are unlikely to be identified, the use of τ b for differential stress estimate is not recommended in practice.

Common Limitations of Calcite Twin Inversion Techniques for Stress and Ways for Improvement
In the future, the penalization function used in CSIT and CSIT-2 could (should?) be improved in a way similar to those suggested by [58,67], and tested on numerically generated or natural polycrystals. A thoughtful comparison among methods applied to the same experimentally deformed samples for which the stress conditions are controlled would usefully complement the theoretical criticisms.
Parlangeau et al. [60] reported that the main technical limitations of paleostress reconstructions from calcite twins to date, whatever the inversion technique used, are related to potential errors in optical measurements. This limitation should be overcome by twin data acquisition using EBSD [72] (Section 4.3).
On top of that, it is strongly advised to take great care when handling the results provided by any stress inversion technique, and to have a critical look at the results, making the best use of available geological information (e.g., geological setting, relative chronology between meso-and macrostructures) to guide and support the separation of superimposed stress tensors. In fault-slip analysis [110], worried about the uncritical use of stress inversion software that is considered by (too) many people as "black boxes" requiring only a few and poor quality fault-slip data as input and yielding indisputable "bulk" paleostress tensors as output. Lacombe [110] claimed that the popularity of such "easy-to-use" computer-based inversion techniques has caused the oversight of the fundamental assumptions they rely upon. Likewise, we would like to call for caution when using calcite twin-based stress inversion techniques since there is a risk-even larger than for fault-slip data-that people uncritically use such populated techniques as black boxes.

Some Regional Studies Using Analysis of Calcite Twins for Paleostress
Interestingly, the authors who were the most critical with the methods predating theirs [58,66] are those whose methods have never or rarely been successfully applied to any experimentally or naturally deformed samples [58,83,111].
It is of course not of great interest here to summarize all the regional studies that were carried out using calcite twins in the last decade [112][113][114][115][116], some of them in a more or less convincing and rigorous way. We summarize hereafter those carried out using stress inversion techniques only.
Aubourg et al. [117] reported a magnetic fabric study in some folds of the Zagros Simply Folded Belt and emphasized that the direction of early folding layer-parallel shortening as revealed by anisotropy of magnetic susceptibility was consistent with paleostress orientations inferred from inversion of calcite twin data [98]. A similar conclusion was reached by [118] from their study in the Sheep Mountain Anticline (Wyoming, USA).
Hnat et al. [33] used [54] method to calculate compression axes for each twin set and the calcite twin data dynamic analysis of the compression axes [108] to define the mean paleostress orientations to decipher indentation tectonics in the Tennessee salient of the southern Appalachians. Layer-parallel paleostress orientations within the thrust belt and the foreland revealed a fanning pattern of paleostress directions that matches the current indenter geometry of the Blue Ridge to the east. The authors proposed that this pattern was related to the advancing Blue Ridge thrust sheet, with differential thrust displacement creating primary belt curvature instead of buttress-induced secondary curvature during shortening as in the Pennsylvania salient to the north in the Appalachian chain.
Tripathy and Saha [119] carried out paleostress reconstruction from inversion of calcite twin data in the Proterozoic Cuddapah Basin, India. Jaya and Nishikawa [120] combined calcite twin data and fault-slip data using the Multiple Inverse Method [121] to reconstruct paleostress in the East Walanae fault zone, South Sulawesi, Indonesia. One of the most interesting aspects of the latter study was the original application of MIM to twinned and untwinned planes separately. In other words, the authors calculated principal stress orientations from twin data only but also principal stress orientations from untwin data (which therefore correspond to impossible stress axes since they are determined virtually from non-deformation data) and used the latter to refine the determination of the former.
Parlangeau et al. [60] tested the CSIT-2 technique not only on numerically generated polycrystals but also on naturally deformed samples from the Monte Nero anticline (Apennines). Especially, measurements were taken from a sample displaying two cross-cutting bed-perpendicular veins filled with calcite. These veins developed before the onset of folding [103] in response to layer-parallel shortening. In the oldest vein, twin inversion yielded two stress tensors. Tensor 1 corresponds to a strike-slip/compressional stress regime with σ 1 axis trending parallel to the vein and σ 3 axis trending perpendicular to the vein. Tensor 2 corresponds to a strike-slip stress regime with σ 1 and σ 3 axes lying within and perpendicular to the youngest vein, respectively, which is in agreement with the opening of the latter. This means that the oldest vein recorded the stress regimes responsible for both its own opening (tensor 1) and the opening of the youngest vein (tensor 2). The strike-slip stress tensor determined from the youngest vein is similar to tensor 2 from the oldest vein, with σ 1 and σ 3 axes lying within and perpendicular to the youngest vein, respectively, and therefore reflects its opening. As expected, no stress tensor similar to tensor 1 from the oldest vein could be retrieved from the youngest vein. As a result, this study confirms (1) the ability of the stress inversion technique to yield stress results consistent with the kinematics of associated mesostructures, and (2) that a relative chronology between successive twinning stresses can be established through the consideration of the orientation of the computed stress axes (as well as of stress ratio) with respect to vein orientation [9] (see Section 7).
Zheng and Shan [122] applied the [67] technique to derive new constraints on the regional tectonic history and the formation of the tectonically complex Huangling Dome (northern South China). However, to test a new technique, one would have expected it to be applied first to a structurally simple target for which independent stress information was available. Instead, the authors have applied for the first time their technique to natural samples from complex, poorly known structures, and interpreted their stress results in terms of orientation and timing through a loose correlation to far-field stresses transmitted through time from distant plate boundaries. While the theoretical validity of the technique itself is beyond doubt, the regional application looks like chasing two rabbits at once without catching none. Further application to simpler structures should help validate this technique in a more convincing way.
Beaudoin et al. [102] performed paleostress reconstructions in the folded sedimentary layers and underlying faulted basement rocks of the Laramide Rattlesnake Mountain Anticline based on analyses of calcite twins, joints/veins and faults. These independent studies yielded comparable results, which provided a mutual validation of the analysed stress markers. All meso-and microstructures were found to have formed as a result of three main tectonic events: the earlier Sevier thin-skinned contraction, the Laramide thick-skinned contraction forming the fold, and the late Basin and Range extension. The comparison of the Rattlesnake Mountain Anticline with the Sheep Mountain Anticline located on the opposite edge of the Bighorn basin and previously studied by [82,123] allowed the authors to propose a conceptual model for the kinematic evolution of Laramiderelated basement-cored anticlines and to extend the previous paleostress reconstruction at the scale of the Bighorn Basin.
Arboit et al. [104,124] discussed the tectonic evolution in the Khao Khwang fold-thrust belt in central Thailand by combining the analyses of fractures, striated faults, and calcite twins. The authors reconstructed the orientations of the stresses associated with the tectonic evolution of this highly deformed intercontinental orogen. Indeed, calcite twins recorded five successive tectonic events that occurred before, during, and after the Triassic Indosinian Orogeny, which ended with the India-Asia collision. Even at a regional scale and in such a complex region, calcite twin analysis was able to reveal all the tectonic stages, starting with three pre-folding events from late Permian to Early Triassic followed by two post folding events, during the Mid-Late Triassic and the Late Mesozoic-Early Cenozoic ( Figure 9A).
More recently, Kulikowski and Amrouch [105] combined calcite twin analysis with geophysical data (collected from borehole image logs and 3D seismic surveys) to decipher the tectonic evolution of a fully buried intracontinental basin in one of the most prolific hydrocarbon provinces in Australia, the Cooper-Eromanga Basin. This multi-scale analysis helped constrain the paleostress orientations and regimes during the six main tectonic events that followed the Adelaidean rifting, with more detail than in previous studies in this same basin over the last 3 decades: (1) Carboniferous strike-slip stress regime with σ 1 oriented NNW-SSE related to the Alice Springs event; (2) Mid-Permian compressional stress regime with σ 1 oriented NW-SE; (3) Late Permian strike-slip stress regime with σ 1 oriented NE-SW related to the Daralingie event; (4) Late Triassic compressional stress regime with σ 1 oriented E-W related to the Hunter-Bowen event; (5) Late Cretaceous compressional stress regime with σ 1 oriented E-W, and (6) Paleogene compressional event with σ 1 oriented N-S ( Figure 9B). This study was one of the first (if not the first) to show the applicability of integrating geophysical analyses with calcite twin inversion for stress to decipher the tectonic evolution of entirely subsurface provinces at a regional scale.

Dating Calcite Twinning Events
Dating twinning will remain a challenge for years to come. While absolute dating of individual twinning deformation in a grain is to date out of reach, establishing a relative chronology between twin sets at the grain scale is possible but challenging. Cross-cutting or offset relationships between twin lamellae from different sets in a grain can be observed using an electron microscope (e.g., [125]), but such observations cannot be done optically simply by using a U-stage, except for secondary twins formed within earlier thick twins (e.g., Type IV twins of [4]). Only in a few particular cases can such geometrical relationships be observed optically [9]. Even in such favourable cases, assuming that the different twin sets formed under different stress regimes at different times, which may be wrong, the extrapolation of such isolated observations to the aggregate that underwent polyphase deformation to constrain its stress/strain history would require a significant number of consistent observations, which are difficult if not impossible to achieve. A relative chronology between different twin sets related to successive distinct tectonic stresses can however be sometimes established indirectly [9]: considerations of the attitude of computed stress/strain axes (1) with respect to vein orientations-in case of superimposed stress regimes, a tensor consistent with vein opening is likely recorded during (or at the latest stage) of vein opening while other (inconsistent) tensors reflect later, post-opening stress regimes-or (2) with respect to bedding dip (e.g., fold test), or the comparison of the stress record from both matrix and veins sometimes allow to indirectly establish a relative chronology [9].
Recently, LA-ICPMS in situ absolute dating of calcite [126][127][128] has opened a new way to constrain indirectly the age of twinning events. One can bracket the age of twinning events by combining the analysis of calcite twins with absolute dating of the calcite cement, either from matrix, vein or fault rocks from which the measurements were taken. In the case of sequential crystallization of successive calcite cements during the diagenetic history (from eogenesis to telogenesis), the occurrence of successive twinning events through time can be captured by bracketing the age of a given twinning event between the age of the youngest dated cement from which the twinning event was recognized in terms of stress or strain and the oldest dated cement from which it is not. Of course, this is very demanding as this implies the initial recognition of different types of cement using for instance cathodoluminescence (Figure 1E), the measurements of a sufficient number (20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30) of twinned grains from each cement to perform inversion for stress or strain and the absolute dating of the different cements.
Dating calcite cements from both matrix and different sets of cross-cutting veins combined with paleostress-strain analysis would therefore help constrain the age and sequence of twinning events. Note however that it is the entire twinning event at the scale of the sample which would be datable this way, and by no mean the formation of an individual twin lamella in a given grain. The twinning event is likely to have a certain duration, i.e., the duration of the application of the causative tectonic stress allowing for the cumulative development of twins within randomly oriented grains hence more or less favourably oriented to twin. To the best of our knowledge, no other study to date has gone beyond such 'indirect' dating of twinning strain.
Two groups of studies have reported results that paved the way to the dating of twinning strain. Nuriel et al. [129] have constrained the onset and evolution of the Dead Sea transform based on in situ U-Pb dating and strain analyses of fault-related twinned calcites. The authors carried out direct dating of syn-faulting calcites from different inactive fault strands of the Dead Sea transform. The dominant horizontal shortening direction recorded in the dated twinned calcites marks the onset of left-lateral motion along the evolving plate boundary so the authors derive the age of onset of the transform motion between 20.8 and 18.5 Ma within an~10-km-wide distributed deformation zone in southern Israel, with subsequent northward propagation (17.1-12.7 Ma) and the establishment of a well-developed >500-km-long plate-bounding fault. Craddock et al. [130] later analysed a single fault-related calcite vein from the Gishron fault zone. U-Pb ages of the nine calcite fillings of this vein range from 20.37 Ma to~13 Ma, with the oldest calcite filling in the middle, younging outward with shearing between the oldest eight zones. All nine calcite fillings have unique mechanical twinning strain records. Despite some variations across the vein, the overall tectonic shortening is subhorizontal and sub-parallel to the Gishron fault (∼N-S) and the Dead Sea transform plate boundary. These studies rely on the assumption that twinning occurred either during or shortly after calcite crystallized, in the absence of later recrystallization due to fluid circulations (e.g., [131]).
The second case study for which calcite twinning and U-Pb dating were carried out on the same calcite material is the Sheep Mountain Anticline (Wyoming, USA). This anticline is a thrust-related, basement-cored NW-SE striking fold that developed in the eastern Bighorn Basin during the late Cretaceous-Paleogene Laramide contraction. Three main bed-perpendicular vein sets were dated using U-Pb calcite geochronology [127]. WNW-ESE oriented vein of the first set formed by 81-72 Ma, before the Laramide contraction. The two other sets are related to the Laramide contraction: the older set consists of NE-SW oriented veins dated to 72-50 Ma, formed during a NE-SW directed layer-parallel shortening. The younger set (50-35 Ma) comprises veins striking NW-SE and developed in response to outer-arc extension at the fold hinge. Amrouch et al. [82] carried out an extensive calcite twin analysis from the cements of these vein sets and the Paleozoic rock matrix. They derived stress and strain tensors with principal axes consistent with the pre-folding WNW-ESE oriented (Sevier?) contraction, the early-folding and late-folding NE-SW Laramide contraction as well as local NE-SW extension at fold hinge. Assuming that (1) precipitation of the dated cements was nearly synchronous with vein opening, and (2) the opening of the successive vein sets was closely linked to the successive paleostresses/strain reconstructed from calcite twins, one can safely propose that the absolute time range of each twinning event in both the matrix and the vein sets is similar to the age of the mesostructures the kinematics of which is consistent with the twinning event, and reflects the time range of the related tectonic event.
Interestingly, the preservation of a calcite twin strain/stress record from dated vein cements consistent with the kinematics of opening of the successive veins, in turn, precludes any later fluid-related dissolution and recrystallization which would have led to younger ages unrelated to twinning deformation and vein opening (e.g., [131]).

Use of Calcite Twins as a Strain Gauge
Twinning strain in a rock can be quantified by several factors, including twin volume fraction, which is equal to the product of the number of twins per unit volume of rock and the average volume of all the twins [6]. Using the correlation between the product of twin density and twin volume, Groshong [53,75] and Groshong et al. [6] developed the Calcite Strain Gauge Technique (CSGT hereinafter) to measure calcite twin strain naturally deformed rocks (see also [132]). The principle of the CSGT will not be recalled herein, but the history of its discovery has been recently revived by [76] (this special issue).
Application of the CSGT requires the following assumptions to be valid: (1) low temperatures (dominance of Type I and Type II twins), (2) random C-axis orientations of calcite, (3) homogenous strain, (4) coaxial deformation, (5) volume constancy and (6) low strain (<15%) (e.g., [85]). If these conditions are not fully met, the calculated strain tensor could be biased or later modified, making strain interpretation tricky and possibly unreliable.
In the last decade, most applications of CSGT were conducted by Craddock and co-workers. The authors tackled several problems using this technique, from the formation and evolution of individual structures (e.g., Early Eocene South Fork Slide, Wyoming, USA [133]) to the tectonic evolution at the regional scale (Sinian system, China, [134]; Appalachian Piedmont Province, [135]; Cretan detachment, [136]; Caribbean plate, [137]). Some investigations were also carried out at the continental scale [138][139][140]. For instance, Craddock et al. [138] reconstructed shortening strain patterns in the Paleozoic Appalachian-Ouachita-Marathon foreland. They concluded that calcite twin strain patterns in the North American midcontinent are complex, with likely tectonic stress transmission originated from different directions at different times as a result of changes in orogenic convergence. Craddock et al. [140] conducted the study of calcite strain patterns across the Proterozoic Penokean orogen, USA and Canada. Finite and calcite twinning strains along (~1500 km) and across (~200 km) the Penokean belt reveal layer-parallel shortening axes that are parallel across the belt, or parallel to the tectonic transport direction. Vertical strain overprints were found in a zone 50-60 km from the Penokean margin and are consistent with strains associated with the stacking of nappes. Craddock et al. [140] concluded that the Penokean orogeny was collisional in type, like the Appalachians, with no rotations of thrust sheets during shortening.
Paulsen et al. [141] are among the few who reported calcite strain analyses from a rift setting. The authors analysed twinned calcite from veins hosted by Neogene sedimentary and volcanic rocks recovered by drill core samples from the Terror Rift system in the southern Ross Sea, Antarctica. Strain analyses revealed mostly homogeneous coaxial strain characterized predominantly by subvertical shortening which is attributed to twinning at relatively shallow depths in an Andersonian normal faulting stress regime induced by sedimentary and ice sheet loading of the stratigraphic sequence. This subvertical shortening is associated with an average extension direction that is subparallel to other indicators of Neogene extension. The overall paucity of a non-coaxial layer-parallel shortening signal leads the authors to favour a tectonic model where the extensional strain reflects Neogene orthogonal rifting rather than right-lateral transtensional rifting.
One of the most recent studies based on CSGT was conducted by [85]. The authors applied CSGT to constrain calcite twinning strains from syn-faulting calcite gouges from a variety of small-offset fault types (strike-slip, normal and thrust faults) and geological settings. Under the assumption that strain and stress were co-coaxial, the authors derive local strain/stress orientation in the vicinity of, and inside the fault zones. For the strikeslip faults, the twinning shortening strain axes are horizontal and at an angle of 0 • -60 • to the respective fault plane (dextral or sinistral) although in the majority of cases the shortening axis is parallel to the fault strike. For the dip-slip normal faults studied, the maximum principal stress/shortening strain axes (σ 1 /ε 1 ) are not oriented at 45 • to the fault plane but are sub-horizontal and either strike-parallel or strike-normal. Thrust faults preserve shortening strain axes parallel to the dip-slip kinematic direction, within the fault plane and not at 45 • to the principal plane. The study concludes that most of the fault stress-strain field results contradict the stress distribution associated with these types of faults as predicted by Coulomb-Anderson theory and that as for most faults σ 1 lies in the plane of the fault, there is no shear stress along any of the faults studied when they are active. One possible pitfall of such study would be that although the authors focused on small-offset faults likely associated with a low amount of deformation, one can wonder whether the assumption of coaxiality between stress and strain is sound. It must be kept in mind that in the vicinity of fault zones and highly deformed areas in which non-coaxial deformation prevails, axes of infinitesimal strain or strain rates or stress are not parallel to the axes of the finite strain. So although the interest of the study lies in the fact it is the first direct attempt at constraining directly the stress/strain field within a variety of small fault zones of known kinematics, further studies are needed to confirm the results.
It must be kept in mind that since deformation is usually partitioned and accommodated by several mechanisms (see Section 1), the CSGT obviously only provides information on the amount of calcite twinning strain, which is usually a minor part of the bulk strain accommodated by the rock.

Use of Calcite Twins as a Stress Gauge
The principle of the oldest and most widely used paleopiezometer [87] is that in a sample without any preferred crystallographic orientation, the relative percentages of grains twinned on 0, 1, 2 or 3 twin plane(s) depend on the applied (σ 1 − σ 3 ) value. Since this relationship has been experimentally calibrated, knowing these relative percentages in a sample, and under the hypothesis of a constant CRSS for twinning (10 MPa), the order of magnitude of (σ 1 − σ 3 ) can be estimated. This method does not take into account the grain size dependence of twinning.
The paleopiezometric approach of [26] is based on the sensitivity of the twinning incidence, twin volume fraction and twin density to differential stress; in turn, estimates of these parameters are used to infer differential stress magnitudes. The first two criteria are found to be largely dependent on grain size as twinning is easier in large grains. The twin density paleopiezometer (number of twins per mm) is nearly grain-size independent and is, therefore, easier to use ( Figure 10A). Compared to the [87] technique, this paleopiezometer does not require any assumption about the existence and value of the CRSS.
The above 'classical' paleopiezometers do not check the mutual compatibility of measured twin systems and do not allow one to relate the differential stress to a given stress event since principal stress orientations are not determined and the mutual consistency of twin sets from which differential stress values are derived is not checked [3]. One can therefore wonder about the significance of 'bulk' differential stress values obtained in the case of polyphase tectonic evolution since the finite twinning strain results from the cumulative effects of non-coeval, differently oriented superimposed stress/strain events.

Low-T • C Paleopiezometers Based on Calcite Twin Density
The use of the [26] twin density criterion provides reasonable results of (σ 1 − σ 3 ) when applied at high temperature, but leads to overestimates of differential stresses when applied at a low-temperature twinning deformation (see discussion in [31]). As a result, the [26] paleopiezometer best applies to twinning deformation at high temperatures (between 200 and 800 • C) and large strains (7-30%).
Sakaguchi et al. [142] performed a series of triaxial compression tests at 200 • C and under 115 MPa confining pressure of medium-grained turbidite sandstones containing scattered calcite particles (0.5 vol. %, average size 130 µm) under the elastic regime. The twins formed are mostly thin twins (thickness < 1 µm) and the mean twin density (estimated by counting the number of twins in one grain, normalized to a unit length of 1 mm in the direction perpendicular to the measured twin set, with uncertainties related to the true dip of the considered twin set) increases with differential stress and shows a good correlation with the maximum elastic stress applied. The best-fitting relation between twin density and stress is (σ 1 − σ 3 )(MPa) = 21.86 × (twin density) 0.57 ( Figure 10A). The further discrete element-based numerical simulation of uniaxial compression tests on a virtual sandstone performed by [142] shows that despite the complex distribution of inter-particle stress at each stress level, the statistical mean stress corresponds to the applied external stress load. This result confirms that the statistically averaged calcite twin density can be regarded as a good stress indicator in real sandstones, with possible extrapolation to other materials, such as concrete with included sparse calcite grains.
Rybacki et al. [23,38] reported a series of low-strain triaxial compression and highstrain torsion experiments over an extended range of temperature (20 to 350 • C) performed on marble and limestone samples. These experiments aimed at examining the influence of stress, temperature, and strain on the evolution of twin density, the percentage of grains with 1, 2 or 3 twin sets, and the twin width and at comparing the temperatures and stresses in the experiments with those inferred from twin morphology and density, respectively. Considering differential stress below 250 MPa as realistic geological stress/T • C conditions, and a real accuracy of stress estimates of ± 50%, the linear twin density was found to increase as the rock hardens with strain and is approximately related to the peak differential stress, (σ 1 − σ 3 ), by the relation: (σ 1 − σ 3 )(MPa) = (19.5 ± 9.8) × (twin density) 1/2 ( Figure 10A). The piezometer should be independent of grain size, but [23] recall that there exists a necessary threshold stress to activate twinning in the fine-grained limestone, which hardly twin at low stresses. This new paleopiezometer shares the same limitations that the previous ones in that it does not check the mutual compatibility of measured twin systems and does not allow to relate the differential stress to a given stress event since principal stress orientations are not determined.  [23]. For [23], the green line shows the regression fit of all data while the red line shows a fit of data below~250 MPa differential stress. (B) Compilation of differential stress estimates from orogenic settings (strike-slip or compressional stress regimes) determined from calcite twinning paleopiezometry and stylolite roughness paleopiezometry as a function of paleodepth of deformation. Modified and implemented after [91].

Calcite Twinning Paleopiezometry Based on Both Residual Stress and Plastic Strain Determination
Chen et al. [143] extensively studied deformation twinning and residual stress of deformed calcite in a coarse-grained metamorphic marble from the Bergell Alps using synchrotron polychromatic X-ray microdiffraction. The authors especially focused on the strain/stress state associated with intersecting and stopping (abutting) mechanical twins. Residual stresses were derived from the strain tensor that is itself derived from Laue diffraction patterns. In the crossing twins, equivalent strain/stress is higher and wider distributed near the twin planes. The study of the stress on stopping twins shows that the shear stress is almost completely released on the twin planes close to twin boundaries, whereas normal stress is concentrated close to twin planes. The average lattice strains, ranging from several hundred to over one thousand microstrains, suggest 60-120 MPa residual stresses, which may indicate that local stresses during twin formation must have exceeded 100-200 MPa during tectonic deformation at metamorphic conditions. This microfocus Laue diffraction approach, therefore, reveals to have the potential to be used as a paleopiezometer based on both residual stress and plastic strain determination.

Determination of Differential Stress Magnitudes Using Calcite Twin Inversion Techniques
The constraints brought by twinned and untwinned e-planes on the determination of the stress tensor have been evaluated [61]. This theoretical study confirmed that both twinned and untwinned e-plane data are necessary to determine differential stresses, but showed that constraints of untwinned e-planes are much looser than those of twinned eplanes. Yamaji [61] further stated that twinned and untwinned e-plane data lose resolution in differential stresses (σ 1 − σ 3 ) if twin lamellae are formed under stress conditions with (σ 1 − σ 3 ) greater than 50-100 MPa, i.e.,~5-10 τ a (assuming a τ a value of 10 MPa).
The ability of the CSIT-2 to detect, separate and determine stress tensors from monophase and polyphase twin datasets including measurements errors or various grain sizes has been demonstrated by testing it on numerically generated twin datasets. In monophase datasets with a one-grain size class and no bias, the loss of information on the determination of differential stress, when the applied stress is >50 MPa as highlighted by [61], is not observed. However, this effect is observed with monophase datasets with two different grain size classes or with polyphase datasets. In agreement with [61], the maximum misfits on the differential stress value are observed at high applied differential stress (75 MPa) and may reach ±37%. The usual maximum methodological uncertainty for the calculation of differential stresses is however expected to be about ±25-30% [60]. This is in agreement with the statement that the misclassification of twin planes as untwinned planes (biased data) gives distortive effect to the estimation of (σ 1 − σ 3 ) using the discrimination of twinned/untwinned e-planes when twins were formed at (σ 1 − σ 3 ) greater thañ 10 τ a [61].
Clearly, because of the dependence of the CRSS on grain size and strain hardening, the CRSS value is not constant [9,28] and the value retained during paleopiezometric studies should be properly adapted. While rough estimates of twinning strain are routinely feasible by considering the density of twins and their thickness, taking into account the grain size is a much difficult task. This is because in a given naturally deformed aggregate, the mean grain size is often variable and the determination of the sizes of the analysed grains can routinely be made on 2D surfaces (thin sections) only. As reported earlier, the 2D appearance of a grain can be misleading since an apparent small grain could correspond to a small section of a larger grain.
In some recent calcite-twin based paleostress analyses, the grain size has been tentatively taken into account. Inversion for stress orientations and stress ratio was carried out separately of the defined (apparent) grain size classes to ensure all classes return the same reduced stress tensors. In a second step, distinctive values of the CRSS τ a as deduced from [13] curves ( Figure 4C) are put into the inversion process for each class to derive differential stress magnitudes, which obviously should be similar for a given state of stress across the whole sample. This approach has successfully been conducted by [60] on samples from the Apennines.
The main limitations of differential stress estimates from calcite twins to date, whatever the inversion technique used, are related to our insufficient knowledge of the variations of the CRSS value for twinning with grain size and grain size distribution in the deformed aggregates. This limitation would be possibly circumvented by defining several grain size classes and treating them separately, keeping in mind that defining true grain size from thin sections remains an unsolved issue. Dividing the whole dataset into sub-classes of homogeneous grain size and performing separately inversion onto homogeneous grain size classes for each of which a constant CRSS value can be adopted provided the number of twin data is sufficient to secure a reliable stress tensor computation) is the best way to reliably define differential stresses [60]. However, taking into account routinely the influence of the grain size distribution (e.g., small/large grains surrounded by large/small grains), which could well impact stress transmission hence stress homogeneity across the aggregate, is unfortunately out of reach to date.
Finally, characterizing more precisely the way the CRSS varies with changing grain size remains a challenge for forthcoming studies. One way to address this problem would be first to synthetize calcite aggregates of homogeneously sized grains and devoid of any pre-existing twins, for instance by compressing and heating calcite powders of controlled grain size. The second step would be to carry out (and to monitor using EBSD and highspeed cameras) deformation experiments on such artificial samples for various grain sizes and under controlled stress conditions to track the onset of twinning at the scale of the grains and of the aggregate. The third step would be to carry out inversion of calcite twin data for differential stress from these samples using the CRSS determined during step 2 and to discuss the reliability of the paleopiezometric results against the applied differential stress.

Case Studies
Rybacki et al. [38] investigated the microstructures of four core samples from the San Andreas Fault Observatory at Depth (SAFOD) using optical and transmission electron microscopy. Sample microstructures indicate intense shearing and dissolution-precipitation as main deformation processes. The samples also contain abundant veins filled with calcite that exhibits a different degree of deformation with evidence for twinning and crystal plasticity. Dislocation density and twin density were used as paleo-piezometers. While possibly reflecting chronologically different maximum paleostress states and/or grain-scale stress perturbations (e.g., co-seismic transient high stress), mean differential stress values from dislocation density 68 ± 46 MPa may represent a lower bound while mean differential stress values from twin density 168 ± 60 MPa represent an upper bound for the paleo-state of stress. The lower stress bound is in line with stress estimates from borehole breakout measurements performed in the pilot hole. Under the assumption of hydrostatic pore pressure and a low intermediate principal stress close to the overburden stress (transpressional stress regime), these data constrain frictional sliding of the San Andreas Fault at the SAFOD site to friction coefficients between 0.24 and 0.31, which supports the hypothesis of a weak San Andreas fault. These differential stress values are consistent with the first order with the compilation by [91].
Another application of [23] low T • C calcite twin density paleopiezometer is provided by [144]. The authors carried out an analysis of microstructures in hydrothermal veins from the sedimentary cover and the igneous basement recovered from Hole U1414A, IODP Expedition 344, to constrain the deformation conditions and tectonic history of the subducting Cocos plate offshore Costa Rica. The application of [23,26] paleopiezometers revealed differential stress magnitudes of 49 ± 11 MPa and 69 ± 30 MPa, respectively, for a maximum depth~470 m (bsf). Deformation temperatures are restricted to the range 170-220 • C, due to the characteristics of the existing twins and the lack of high-temperature intracrystalline deformation mechanisms. These differential stresses likely reflect intraplate stresses due to the bending of the plate at the trench [144]. Although determined from a different tectonic setting and using a different method ( [23] vs. [87]), these differential stress results well compare to differential stresses (~48 MPa) associated with the twinning of calcite veins and amygdule fillings within Icelandic basalts and thought to represent ridge-push tectonic conditions [145]. However, the values are slightly higher than the 26 MPa intraplate horizontal differential stresses determined from calcite from amygdules and veins in the lower basaltic portions of DSDP Hole 433 (depth~550 m) in Suiko seamount of the Pacific plate, possibly related with hotspot plumes and/or transmitted from distant plate boundaries [146]. It must be emphasized that calcite twinning is to date the only paleopiezometer that has provided information on differential stress magnitudes in different domains of oceanic plates (ridge-Iceland, intraplate-Pacific), and near plate boundary-offshore Costa Rica).
As a follow-up to the pioneering paper by [82,102] provided integrated pictures of stress magnitudes associated with fold development using calcite twinning paleopiezometry. The authors determined differential stress magnitudes related to the Laramide compressional event at the Rattlesnake Mountain anticline (western Bighorn basin, USA) and compared them to those reconstructed at the Sheep Mountain anticline in the eastern part of the basin by [82], after normalizing them to a similar deformation depth following [3]. Interestingly, the trend of evolution of differential stress (σ 1 − σ h ) (σ h being either σ 3 or σ 2 for a strike-slip or compressional regime, respectively) associated with early-folding layer-parallel shortening and late-stage fold tightening in the backlimb and the forelimb of each fold structure revealed to be similar for both folds. Note however that stress perturbation above the basement thrust as identified at Sheep Mountain anticline during layer-parallel shortening [82] could not be retrieved at Rattlesnake Mountain anticline presumably because the formations located in a similar structural position than in Sheep Mountain Anticline have been eroded away. As a result, Laramide differential stress magnitudes were seemingly not attenuated across the Bighorn Basin and is likely controlled by stress transmission throughout the basement during thick-skinned tectonics rather than by the distance to the orogenic front as in thin-skinned tectonics settings. This conclusion is in agreement with that derived from the later combination of calcite twinning and stylolite roughness paleopiezometers at Sheep Mountain anticline [147] (see Section 10.2).
In addition to the study of contemporary fault stability and/or reactivation under the modern state of stress using geomechanical analysis [148], calcite twin-based paleostress determination may provide useful insights on fault and fracture evolution through time in regions with poor seismic imaging. Kulikowski and Amrouch [149] performed a 4D quantification of stress in the Cooper-Eromanga basin as well as estimates of fault reactivation, dilation and slip tendency ( Figure 11) using an approach combining inversion of calcite twins for stress orientations and magnitudes with rock mechanics. Kulikowski and Amrouch [149] were indeed able to estimate fault reactivation potential at key stages of the petroleum system development. This kind of approach has the potential to provide valuable information for subsurface resources exploration, by upgrading the existing static models of fluid migration into dynamic ones, taking into account the 4D structural permeability behaviour. Roure et al. [150] reviewed the use of paleo-thermo-barometers (including calcite twinning paleopiezometry) and coupled thermal, fluid flow and pore-fluid pressure modelling for hydrocarbon and reservoir prediction in fold-and-thrust belts.
Last, original applications of calcite twinning paleopiezometry to planetary sciences have been provided by [151,152]. Lindgren et al. [151] used calcite twinning paleopiezometry to constrain the pressure threshold of impacts on extraterrestrial bodies, while [152] study aimed at constraining subsurface deformation of experimental hypervelocity impacts in quartzite and marble targets. Calcite e-twins were estimated as starting to form at shock pressures of ∼110 to 480 MPa as inferred from modelling [151], which is at least a factor of ten higher than the 10 MPa that is considered to be required to produce calcite twins in low strain rate terrestrial settings. Whereas these pressure estimates are consistent with calcite twinning in carbonaceous chondrites being a result of impact gardening in shallow levels of their asteroidal parent bodies, such estimates should be considered with caution considering the poor knowledge we have on the actual deformation behaviour of calcite at such high pressure/high-velocity conditions. at shock pressures of ∼110 to 480 MPa as inferred from modelling [151], which is at least a factor of ten higher than the 10 MPa that is considered to be required to produce calcite twins in low strain rate terrestrial settings. Whereas these pressure estimates are consistent with calcite twinning in carbonaceous chondrites being a result of impact gardening in shallow levels of their asteroidal parent bodies, such estimates should be considered with caution considering the poor knowledge we have on the actual deformation behaviour of calcite at such high pressure/high-velocity conditions. Figure 11. Paleostress parameters evolution through time in the Cooper-Eromanga Basin as deduced from inversion of calcite twin data. The reactivation potential of the four most common fault sets measured in the Cooper-Eromanga Basin is estimated by calculating the dilation tendency and the slip tendency. Modified after [149].

Combining Calcite Twinning Paleopiezometry with Fractures and Rock Mechanics Data
For a given tectonic event, the reconstruction of the effective principal stress magnitudes may be achieved by combining paleopiezometry from calcite twins with analysis of Figure 11. Paleostress parameters evolution through time in the Cooper-Eromanga Basin as deduced from inversion of calcite twin data. The reactivation potential of the four most common fault sets measured in the Cooper-Eromanga Basin is estimated by calculating the dilation tendency and the slip tendency. Modified after [149].

Combining Calcite Twinning Paleopiezometry with Fractures and Rock Mechanics Data
For a given tectonic event, the reconstruction of the effective principal stress magnitudes may be achieved by combining paleopiezometry from calcite twins with analysis of coeval fractures (faults and/or veins) and mechanical properties of rocks as derived from mechanical tests [3,8,97,153]. The basic principle is to find the values of σ 1 , σ 2 and σ 3 required for consistency between newly formed faulting, frictional sliding along pre-existing fault planes and calcite twinning. Knowing the differential stress value (σ 1 − σ 3 ) from calcite twinning paleopiezometry (i.e., the diameter of the Mohr circle), the complete stress tensor can be determined by fitting the Mohr circle to either the failure envelope (in case of newly formed faults) or to the Byerlee's friction line (for reactivated inherited faults), or both (see also [91,154]).
A step forward consists of comparing the magnitude of the effective vertical principal stress determined with the above approach with the vertical stress magnitude independently determined from overburden at the time of twinning and faulting, assuming a hydrostatic fluid pressure ( Figure 12). In the absence of variation of the vertical stress due to erosion or deposition, the difference between the two values likely reflects fluid overpressure ( Figure 12). This approach applied to the Laramide Sheep Mountain and Rattlesnake Mountain anticlines have provided unprecedented information on the stress history and evolution of pore-fluid pressures during folding and contraction in the Bighorn Basin (Wyoming, USA). Especially, fluid pressure was found (1) to increase up to lithostatic during Laramide layer-parallel shortening, possibly due to either porosity reduction by pressure-solution, poor hydraulic permeability of fracture sets due to their low vertical persistence or to their fast healing, strong increase in horizontal stress magnitude or input of exotic fluids into the reservoir in response to a large-scale fluid migration, then (2) to drop down back to~hydrostatic during folding due to development of curvature-related fractures that enhanced the hydraulic permeability of the reservoir [123,155]. This break of fluid compartmentalization within the formations of interest is consistent with the geochemical signatures of syn-folding vein cements, which suggest a vertical migration of deeper radiogenic hot fluids within the sedimentary cover during folding [156].
coeval fractures (faults and/or veins) and mechanical properties of rocks as derived from mechanical tests [3,8,97,153]. The basic principle is to find the values of 1, 2 and 3 required for consistency between newly formed faulting, frictional sliding along pre-existing fault planes and calcite twinning. Knowing the differential stress value (σ1 − σ3) from calcite twinning paleopiezometry (i.e., the diameter of the Mohr circle), the complete stress tensor can be determined by fitting the Mohr circle to either the failure envelope (in case of newly formed faults) or to the Byerlee's friction line (for reactivated inherited faults), or both (see also [91,154]).
A step forward consists of comparing the magnitude of the effective vertical principal stress determined with the above approach with the vertical stress magnitude independently determined from overburden at the time of twinning and faulting, assuming a hydrostatic fluid pressure ( Figure 12). In the absence of variation of the vertical stress due to erosion or deposition, the difference between the two values likely reflects fluid overpressure ( Figure 12). This approach applied to the Laramide Sheep Mountain and Rattlesnake Mountain anticlines have provided unprecedented information on the stress history and evolution of pore-fluid pressures during folding and contraction in the Bighorn Basin (Wyoming, USA). Especially, fluid pressure was found (1) to increase up to lithostatic during Laramide layer-parallel shortening, possibly due to either porosity reduction by pressure-solution, poor hydraulic permeability of fracture sets due to their low vertical persistence or to their fast healing, strong increase in horizontal stress magnitude or input of exotic fluids into the reservoir in response to a large-scale fluid migration, then (2) to drop down back to ~hydrostatic during folding due to development of curvature-related fractures that enhanced the hydraulic permeability of the reservoir [123,155]. This break of fluid compartmentalization within the formations of interest is consistent with the geochemical signatures of syn-folding vein cements, which suggest a vertical migration of deeper radiogenic hot fluids within the sedimentary cover during folding [156]. Figure 12. Principle of the determination of the absolute principal stress magnitudes using inversion of calcite twin data combined with fracture analysis and rock mechanics data [153] and of assessment of the fluid overpressure [123]. Vth is the vertical stress assuming a hydrostatic fluid pressure. Veff is the actual vertical stress. V is the difference between Vth and Veff. Modified after [155]. Figure 12. Principle of the determination of the absolute principal stress magnitudes using inversion of calcite twin data combined with fracture analysis and rock mechanics data [153] and of assessment of the fluid overpressure [123]. σ Vth is the vertical stress assuming a hydrostatic fluid pressure. σ Veff is the actual vertical stress. ∆σ V is the difference between σ Vth and σ Veff . Modified after [155].

Combining Calcite Twinning Paleopiezometry with Stylolite Roughness Paleopiezometry
The main successful and promising progress has been to combine calcite twinning paleopiezometry with stylolite roughness paleopiezometry [91,103,147]. Figure 13A,B illustrates the principle of the combination of calcite twinning paleopiezometry and stylolite roughness paleopiezometry and its possible derivatives.

Combining Calcite Twinning Paleopiezometry with Stylolite Roughness Paleopiezometry
The main successful and promising progress has been to combine calcite twinning paleopiezometry with stylolite roughness paleopiezometry [91,103,147]. Figure 13A,B illustrates the principle of the combination of calcite twinning paleopiezometry and stylolite roughness paleopiezometry and its possible derivatives. Figure 13. (A). General principle of the combination of calcite twinning and stylolite roughness paleopiezometers with burial history of strata to derive the timing of the successive deformation stages in contractional domains (e.g., fold-andthrust belts). (B) Typical tectonic association and relative chronology between sedimentary stylolites, tectonic veins and tectonic stylolites as observed in the field. Principle of stylolite roughness paleopiezometry applied to sedimentary stylolites (C) and to tectonic stylolites (D) and related outputs. (E) Principle of the combination of calcite twinning paleopiezometry and stylolite roughness paleopiezometry applied to tectonic stylolites to constrain the time at which layer-parallel shortening preceding folding stopped. For more details, see the text.
The principle of the stylolite roughness paleopiezometry is summarized hereinafter. Stylolites are rough surfaces of localized dissolution developed in the rock under applied pressure. Stylolite teeth are parallel to the maximum principal stress σ1, related either to compaction (sedimentary stylolites, vertical σ1) or to tectonics (tectonic stylolites, horizontal σ1). As the name hints, stylolite roughness paleopiezometry builds on signal properties of the stylolite roughness along 2D tracks. The roughness is defined as the difference in height between two points of a track, separated by a set given length, defining the socalled scale of observation. Empirical observations conducted by [157] showed that the analysis of the roughness returns two regimes of self-affine properties, according to the scale of observation. In case the roughness is considered over a large-scale, i.e., when the along-track length separating the two considered points is typically above 1 mm, the relationship between the transform of the roughness and the considered length or spatial frequency is such that the roughness coefficient (so-called Hurst exponent) equals 0.5, a value typical of a control by the elastic energy. When the roughness is considered with a smaller scale, with a length between two points below typically 1 mm, the related roughness coefficient is ~1.1, suggesting a control by the surface energy [157]. Such finding led [158] to propose an analytical relationship linking the applied differential stress, the mean stress, the elastic properties of the rock (Young modulus and Poisson coefficient) and the surface energy at the solid-fluid interface to the characteristic length (referred to as the cross-over length, Lc) at which roughness analysis switches from a Hurst coefficient related to elastic energy to a Hurst coefficient related to surface energy ( Figure 13C). Since then, the Stylolite Roughness Inversion Technique (SRIT) has successfully been applied to The principle of the stylolite roughness paleopiezometry is summarized hereinafter. Stylolites are rough surfaces of localized dissolution developed in the rock under applied pressure. Stylolite teeth are parallel to the maximum principal stress σ 1 , related either to compaction (sedimentary stylolites, vertical σ 1 ) or to tectonics (tectonic stylolites, horizontal σ 1 ). As the name hints, stylolite roughness paleopiezometry builds on signal properties of the stylolite roughness along 2D tracks. The roughness is defined as the difference in height between two points of a track, separated by a set given length, defining the so-called scale of observation. Empirical observations conducted by [157] showed that the analysis of the roughness returns two regimes of self-affine properties, according to the scale of observation. In case the roughness is considered over a large-scale, i.e., when the alongtrack length separating the two considered points is typically above 1 mm, the relationship between the transform of the roughness and the considered length or spatial frequency is such that the roughness coefficient (so-called Hurst exponent) equals 0.5, a value typical of a control by the elastic energy. When the roughness is considered with a smaller scale, with a length between two points below typically 1 mm, the related roughness coefficient is~1.1, suggesting a control by the surface energy [157]. Such finding led [158] to propose an analytical relationship linking the applied differential stress, the mean stress, the elastic properties of the rock (Young modulus and Poisson coefficient) and the surface energy at the solid-fluid interface to the characteristic length (referred to as the cross-over length, Lc) at which roughness analysis switches from a Hurst coefficient related to elastic energy to a Hurst coefficient related to surface energy ( Figure 13C). Since then, the Stylolite Roughness Inversion Technique (SRIT) has successfully been applied to sedimentary stylolites from various geological settings, from simple ones [159][160][161] to more complex ones [103,[162][163][164][165]. Assuming that sedimentary stylolites formed under a vertical maximum principal stress σ 1 with the horizontal stress being isotropic (σ 2 = σ 3 ), applying SRIT to sedimentary stylolites provides a reliable way to access the magnitude of the past vertical stress, hence the burial depth, at the moment the stylolite stop dissolving, free from any considerations upon fluid pressure or geothermal gradient [91,160] (Figure 13B,C). When applied to a population of sedimentary stylolites, SRIT yields the range of depths in which compaction-induced pressure solution was active under a vertical σ 1 .
On the other hand, fewer attempts at applying SRIT to tectonic stylolites have been published [103,147,166]. Yet the technique seems perfectly valid beyond some methodological issues related to the fact that the state of stress is by essence triaxial, with σ 1 being generally horizontal and normal to the stylolite at the time of its formation, leading to a stress anisotropy within the stylolite plane (σ 2 = σ 3 ) that translates into a varying value of Lc with respect to the in-plane orientation of σ 2 and σ 3 ( Figure 13D). The relationship between Lc and the angle θ is a periodic function, with minimum and maximum Lc separated by 90 • . Methodologically, when not studying open stylolites 3-D surfaces [166], using tectonic stylolites requires to conduct roughness inversion on 2D scans of at least three surfaces normal to the stylolite corresponding to different in-plane angles θ. This procedure returns three Lc corresponding to the three corresponding angles between the cuts and the vertical direction, allowing for the reconstruction of the Lc variations in 3D. The minimum and the maximum Lc correspond to (σ 1 − σ 3 ) and (σ 1 − σ 2 ), respectively ( Figure 13D). If σ v at the time of the tectonic stylolites development is known independently, SRIT yields principal stress orientations and absolute magnitudes of horizontal principal stresses.
If applied separately, stylolite roughness paleopiezometry and calcite twinning paleopiezometry may provide a record of the progressive stress loading during layer-parallel shortening of flat-lying strata leading ultimately to faulting and/or folding. One can for example expect that during stress loading, pressure solution initiates at low differential stress along planar solubility heterogeneities in rocks and stops rapidly, presumably by clogging around the dissolution planes [167]. Stylolites would therefore dissipate a low amount of stress before saturating. Increasing differential stress would then trigger later calcite twinning, then fracturing under higher differential stress. In this scenario, stylolite development would mostly predate calcite twinning strain, which seems to be confirmed by recent studies [147]. An important derivative would therefore be a better understanding of the timing and significance of the stress record by the distinctive paleopiezometers. However, a difficulty arises which is related to the time-scale at which the various paleopiezometers record stress, i.e., progressive vs. instantaneous stress record. At our present state of knowledge, it seems that the calcite twinning paleopiezometer likely provides differential stress reached following a long-lasting stress loading during a given tectonic event, while stylolite roughness paleopiezometry applied to a population instead yields snapshots of stress at different times of ongoing deformation [91].
One can wonder whether, in addition to recording the state of stress at a given, instantaneous time, tectonic stylolites would also record a very local state of stress only. If this is the case, this would probably limit the interpretation of the derived absolute stress magnitudes to the local stress heterogeneities rather than to regional trends. It would therefore be complicated to use results from SRIT alone to unravel large-scale governing parameters for stress distribution in the upper crust.
The combined application of calcite twinning and stylolite roughness paleopiezometric techniques has provided a partial answer to this question, as well as a mutual validation of the stress results obtained using both paleopiezometers. Beaudoin et al. [103] have reported on the complex paleostress and paleoburial history of the Monte Nero Anticline (Apennines). The analysis of sedimentary and tectonic stylolites, veins/fractures and calcite twins enabled the reconstruction of the principal stress orientations and magnitudes associated with the main stages of fold evolution, i.e., layer-parallel shortening, fold growth and related extension at fold hinge and late-stage fold tightening. The absolute stress magnitudes determined from calcite twinning and stylolite roughness paleopiezometry, as well as their evolutionary trend through time during orogenic contraction, are comparable ( Figure 14A). Interestingly, the authors further documented a local switch of the stress regime from contractional to extensional as a result of erosion and local sediment redistribution or structural burial in the overturned forelimb of the fold, along with switches between reverse and strike-slip stress regimes during layer-parallel shortening [91] (Figure 14A). By allowing for a quantification of stress magnitudes and by documenting an unexpected stress complexity during contractional deformation, this approach goes beyond classical views [168,169] and provides new insights into the stress evolution in fold-and-thrust belts.
Yet the combined application of CSIT and SRIT to the Sheep Mountain-Little Sheep Mountain Anticlines in the Bighorn Basin enabled [147] to quantify and compare the stress records associated with the distinctive tectonic styles of deformation, namely thin-skinned and thick-skinned. The folds formed as basement-cored Laramide structures in the formerly undeformed foreland of the thin-skinned Sevier orogen. Differential stress magnitudes that prevailed during Sevier and Laramide layer-parallel shortening that preceded the onset of large-scale Laramide folding were determined independently from stylolite roughness and calcite twinning paleopiezometry in the same sedimentary formations. Differential stress transmitted from the distant Sevier thin-skinned orogeny into the sedimentary cover of the Bighorn basin was found to be higher than that associated with the early stage of LPS related to the thick-skinned Laramide deformation. This suggests that the Sevier stress recorded in the Bighorn Basin was efficiently transmitted from the distant Sevier thin-skinned orogeny into the stable foreland through a shallow stress guide, i.e., the sedimentary cover ( Figure 14B). In contrast, the Laramide stress was likely transmitted forelandward through a deep crustal stress guide, so most of the stress would have been dissipated at depth while triggering the inversion of inherited basement normal faults and only a fraction of this stress was transmitted upward into the attached cover during the early stage of Laramide contraction [147]. This study, therefore, shows that the tectonic style of an orogen does affect the transmission of early orogenic stress into the stable continental interior, thus challenging the conclusions of earlier studies [170].
Ongoing and future development combining both paleopiezometers mainly rely on the capability of SRIT to constrain the σ v prevailing at the time of deformation. A first promising development is to better constrain the timing at which σ 1 switches from a vertical to a horizontal attitude by combining the SRIT results from a population of sedimentary stylolites and burial models ( Figure 13A, [164,165]). This approach constrains the range of depths at which twinning occurred under a horizontal σ 1 , along with providing indirectly the absolute timing of the contractional event ( Figure 13A). Future development can use the mathematical relationship between the horizontal differential stress (σ 1 − σ h ) determined from tectonic stylolites and the magnitude of σ v to better estimate the depth of deformation. Assuming that despite a possible distinctive stress record (see above) the maximum differential stress (σ 1 − σ h ) provided independently by CSIT and SRIT should be roughly consistent, one can perform an iterative search for the maximum value of σ v that may allow for this consistency ( Figure 13E) and which presumably corresponds to the maximum depth at which layer-parallel shortening was recorded in the rock. When reported on the burial-time curve of the formation from which micro-and meso-structures were collected ( Figure 13A), this depth estimate provides a timing for the end of layerparallel shortening, just before (or at the very onset of) exhumation related to fold growth started. This approach may pave the way to a better appraisal of the burial-stress-time evolution at the reservoir scale and the fluid pressure. Other developments can arise from the study of tectonic, sedimentary-like stylolites related to local or regional extension. Indeed, bedding-parallel stylolites are often documented near normal faults (e.g., [171]), and treating those as tectonic stylolites, i.e., by reconstructing a periodic anisotropy of Lc (see above), could help characterize the local state of stress, including absolute stress magnitudes. To do so, an independent estimate of the extensional differential stress is required, motivating again a combination with calcite twins paleopiezometry.
Geosciences 2021, 11, x FOR PEER REVIEW 41 of 48 above), could help characterize the local state of stress, including absolute stress magnitudes. To do so, an independent estimate of the extensional differential stress is required, motivating again a combination with calcite twins paleopiezometry. Figure 14. (A) Absolute stress magnitude of the vertical stress (v, circles), the minimum horizontal principal stress (h, squares) and the maximum horizontal principal stress (H, crosses) plotted against the timing of corresponding deformation. Data come from the Monte Nero Anticline, Italy after [91,103]. Timing is estimated after [164]. The top plot (cyan symbols) reports the magnitudes reconstructed from CSIT, bottom plot (red symbols) reports the magnitudes reconstructed from SRIT. Note that the magnitudes of both horizontal stresses are unknown during the folding phase of the fold development, yet they are locally inferior to the vertical stress. (B) Differential stress magnitudes (σ₁ − σ₃) determined from inversion of tectonic stylolite roughness (squares) and calcite twins (circles) at Sheep Mountain and Little Sheep Mountain anticlines for the Laramide (green) and Sevier (blue) events. The insert summarizes a conceptual model of stress transmission through shallow and deep stress guides that accounts for the difference in differential stress magnitudes sustained by sedimentary cover rocks during layer-parallel shortening related to thin-skinned Sevier orogeny and thickskinned Laramide orogeny. Modified after [147]. Data come from the Monte Nero Anticline, Italy after [91,103]. Timing is estimated after [164]. The top plot (cyan symbols) reports the magnitudes reconstructed from CSIT, bottom plot (red symbols) reports the magnitudes reconstructed from SRIT. Note that the magnitudes of both horizontal stresses are unknown during the folding phase of the fold development, yet they are locally inferior to the vertical stress. (B) Differential stress magnitudes (σ 1 − σ 3 ) determined from inversion of tectonic stylolite roughness (squares) and calcite twins (circles) at Sheep Mountain and Little Sheep Mountain anticlines for the Laramide (green) and Sevier (blue) events. The insert summarizes a conceptual model of stress transmission through shallow and deep stress guides that accounts for the difference in differential stress magnitudes sustained by sedimentary cover rocks during layer-parallel shortening related to thin-skinned Sevier orogeny and thick-skinned Laramide orogeny. Modified after [147]. Figure 10B compiles differential stress estimates from orogenic settings (horizontal σ 1 , strike-slip or compressional stress regimes) determined from calcite twinning paleopiezometry and stylolite roughness paleopiezometry as a function of paleodepth of deformation. It therefore encompasses only part of the previous release by [91] (to which the reader can refer for the origin of stress estimates) but has been implemented with more recent estimates, especially using SRIT. At the first order, differential stress increases with depth, which supports a long-term frictional behaviour of the upper crust, with higher stress, hence higher crustal strength, in compressional regime than in strike-slip regime. Interestingly, most differential stress data plot along, or close to, the stress-depth curves predicted for a critically stressed crust for a range of friction coefficients and pore pressure ratios [3,8,91]. Outliers, most of those coming from stylolite roughness inversion, may reflect either long-term stress levels kept beyond the frictional yield thanks to stress release, e.g., by pressure-solution creep (upper-tier of the ductile flow regime [172], or more likely transient stress levels recorded as snapshots during stress loading which ultimately leads to calcite twinning and frictional faulting.
Again, it must be kept in mind that the compiled estimates are only those from which the corresponding paleodepth of deformation could be evaluated independently. Assessing the depth at the time of the deformation is a key parameter in paleopiezometry. Differential stress estimates from different sources and different settings cannot be reliably compared if not normalized to the same depth of deformation.

Conclusions
This review reports the significant progress that has been made over the last decade in (1) the understanding of the occurrence of calcite deformation twins, (2) the methodological development of both acquisition and treatment of calcite twins and (3) strain and stress gauge applications. It proposes a critical snapshot of the state of the art, highlighting some new findings partly challenging established concepts, such as paleothermometry based on twin thickness. This overview also points out the promising future research, which should focus on the better appraisal of the meaning and level of stress required for twins to nucleate and grow. Depending on the purpose of calcite twin data collection (e.g., fabric analysis in marbles and strain analysis vs. stress analysis), measurement techniques and procedures could be adapted, keeping in mind the necessary balance between the requirements for a sufficient number of grains and for the accuracy of measurements versus measurement error and the time and cost of data acquisition. Thanks to these recent developments, the semi-automated calcite twin data acquisition is no longer out of reach.
The application of calcite twinning paleopiezometry to tectonic studies and its combination with other paleopiezometric techniques have also been a significant step forward. It helped address important questions on the controls of the level of stress in the shallow crust and how this stress evolves in time and space in various settings and in response to various conditions. Improving our knowledge of the stress and its variations across scales is of prime importance for applied geological purposes such as geological hazards, engineering activities and resource exploration, but also for fundamental geological purposes, such as understanding the mechanical behaviour of geological materials and deciphering various tectonic mechanisms. Calcite twinning paleopiezometry remains to date one of the most precious tools to better understand and quantify the stress evolution in the upper crust. Stress estimates from calcite twins will help feed larger-scale geomechanical models, with important implications for the appraisal of crustal rheology, earthquake processes or fluid distribution in reservoir rocks. We believe that alone must motivate the community to keep applying and developing these techniques.