Evolution of a hydrothermal ore-forming system recorded by sulfide mineral chemistry: a case study from the Plaka Pb–Zn–Ag Deposit, Lavrion, Greece

Laser ablation-inductively coupled plasma-mass spectrometry and electron-probe microanalysis were used to investigate the trace-element contents of sphalerite, chalcopyrite and pyrite from the Plaka Pb–Zn–Ag deposit. Using petrographic observations, the analytical results could be linked to the temporal evolution of the Plaka ore-forming system. Sphalerite chemistry reliably records the temperature and fS2 evolution of the system, with estimated formation temperatures reproducing the microthermometric results from previous fluid-inclusion studies. Chalcopyrite chemistry also shows systematic variations over time, particularly for Cd, Co, Ge, In, Sn and Zn concentrations. Measurable pyrite was only found in association with early high-temperature mineralisation, and no clear trends could therefore be identified. We note, however, that As and Se contents in pyrite are consistent with formation temperatures estimated from co-existing sphalerite. Statistical analysis of the sphalerite data allowed us to identify the dominant geological controls on its trace-element content. The three investigated factors temperature, fS2, and sample location account for > 80% of the observed variance in Mn, Fe, Co, Ga, Ge, In, Sb and Hg concentrations, and > 60% of the observed variance in Cd and Sn concentrations. Only for Cu and Ag concentrations is the explained variance < 50%. A similarly detailed analysis was not possible for chalcopyrite and pyrite. Nevertheless, comparison of the results for all three investigated minerals indicates that there are some systematic variations across the deposit which may be explained by local differences in fluid composition.


Introduction
The trace-element signatures of sulfide minerals are quickly becoming important tools in the study of hydrothermal mineral deposits Fontboté et al. 2017). This is because they can record critical information on both the physical and chemical conditions of ore formation, as demonstrated by several recent studies (Deditius et al. 2014;Frenzel et al. 2016Frenzel et al. , 2020Keith et al. 2018). However, the field is still young, and relatively little work has been done to systematically link changes in sulfide trace-element chemistry to the evolution of individual ore-forming systems (e.g. Sykora et al. 2018;Bauer et al. 2019;Godefroy-Rodriguez et al. 2020). Studies which simultaneously investigate the chemistry of multiple co-existing sulfide minerals are even rarer (e.g. George et al. 2016).
In this work, we used laser ablation-inductively coupled plasma-mass spectrometry (LA-ICP-MS) and electronprobe microanalysis (EPMA) to investigate trace-element concentrations in sphalerite, chalcopyrite and pyrite in Pb-Zn-Ag ores from different parts of the Plaka deposit, Lavrion, Greece. The Plaka deposit is characterised by an evolution from high-to low-temperature sulfide mineralisation styles (Voudouris et al. 2008). This evolution is reasonably well constrained (Voudouris et al. 2008), making it an excellent case study to: (1) examine changes in the traceelement signatures of the investigated minerals with oreforming conditions and (2) test the applicability of existing geochemical tools such as the GGIMFis (= Ga, Ge, In, Mn and Fe in phalerite) geothermometer (Frenzel et al. 2016) at the scale of a single deposit.

District geology
The Lavrion area is located 50 km southeast of Athens ( Fig. 1A) and belongs to the Attic-Cycladic Crystalline Belt, a metamorphic terrain formed between the Late Cretaceous and Late Miocene (Altherr et al., 1982, Katzir et al. 2000Bröcker et al. 2004;Bröcker and Keasling 2006;Scheffer et al. 2016). Two tectonic units dominate the geology of the area. They are separated by a large detachment fault (Fig. 1A). The basal Kamariza Unit is composed of a sequence of Triassic to Early Jurassic metasedimentary rocks. It is subdivided into three members: The Lower Marble, Kamariza Schist and Upper Marble ( Fig. 1A; Marinos and Petrascheck 1956, Photiades and Carras, 2001, Scheffer et al. 2016. The Lavrion Unit overlies the Kamariza Unit and consists of the Lavrion Schist, marbles and a meta-ophiolite (Photiades and Carras, 2001). During peak metamorphism, both units reached blueschist facies conditions (Scheffer et al. 2016). However, they later retrogressed mostly to greenschist facies mineral assemblages (Scheffer et al. 2016).
Several occurrences of granite and granodiorite laccoliths, pipes, dikes and sills are exposed across the Lavrion district (Marinos and Petrascheck 1956;Skarpelis 2007;Papanikolaou and Syskakis 1991;Skarpelis et al. 2008). The most important of these is the Plaka intrusion (Fig. 1B), an I-type granodiorite surrounded by an extensive contact metamorphic aureole of calc-silicate hornfels in the Kamariza schists ( Fig. 1B; Baltatzis 1981). Field relationships indicate that magmatism was broadly synchronous with the post-metamorphic development of the large detachment fault, with available radiometric dates suggesting an absolute age of 9.7-8.1 Ma (Skarpelis et al. 2008;Liati et al. 2009).

Mineralisation
The Lavrion district is famous for its long mining history, starting before 3000 bc (Conophagos 1980;Roald and Webster 2018) and ending in the 1970s (Marinos and Petrascheck 1956). Two major mining centres were located at Plaka and Kamariza, but several other sulfide deposits occur in the district (Fig. 1A). Overall, about 2.3 Mt of Pb, 1.12 Mt of Zn and 7.8 kt of Ag were produced from the district (Conophagos 1980). However, not all Zn originally contained in the ores was extracted, so that a reasonable estimate of the original Zn:Pb ratio in the ores is ~ 1:1.
Mineralisation appears to be linked to Late Miocene magmatism (Marinos andPetrascheck, 1956, Voudouris et al. 2008;Bonsall et al. 2011). This is particularly clear in the Plaka area where the different mineralisation styles developed after hornfels formation in and around the Plaka granodiorite, mostly in the Kamariza Schist and Upper Marble formations (Leleu et al. 1973;Economou et al. 1981;Fig. 1B).
Fluid-inclusion studies indicate a two-stage evolution of the Plaka mineralising system. Early fluids were high-temperature and high-salinity magmatic fluids (~ 360 °C, < 40 wt.% NaCl eq ) that mostly deposited early porphyry-style, skarn and carbonate-replacement mineralisation. These were followed by low-temperature and lowersalinity (< 250 °C, < 10 wt.% NaCl eq ) fluids of mixed origin responsible for the formation of late Ag-rich vein-hosted assemblages (Voudouris et al. 2008;Bonsall et al. 2011;Scheffer et al. 2019).
Sulfide parageneses in all four Pb-Zn-dominated mineralisation styles (skarn-hosted and skarn-free carbonate-replacement, breccia-hosted and vein-type) are similar, starting with early pyrrhotite followed by pyrite, sphalerite, chalcopyrite, galena and sometimes arsenopyrite (Voudouris et al. 2008). This early paragenesis corresponds to the early high-temperature fluids. Vein-type mineralisation additionally contains a later paragenetic stage comprising base-metal sulfides, as well as native arsenic, sulfosalts and various silver minerals (Voudouris et al. 2008). This later paragenesis corresponds to the lower-temperature fluids. Sphalerite occurs across both paragenetic stages and evolved from iron-rich to iron-poor compositions as the system cooled (Voudouris et al. 2008).

Samples
Six samples of carbonate-replacement and vein-style mineralisation were selected for detailed investigation (Table ESM1.1; Fig. 2). The material was collected during an earlier sampling campaign from the surface and underground galleries at Plaka (Voudouris et al. 2008). Sampling locations are shown in Fig. 1B. Polished rounds (ø 25 mm) were prepared from all samples at the National and Kapodistrian University of Athens. The focus during sample selection was on material containing abundant sphalerite associated with other sulfides. Sphalerite is currently the best understood sulfide mineral in terms of the geological controls on its trace-element content (e.g. Lusk and Calder 2004;Frenzel et al. 2016Frenzel et al. , 2020 and was therefore of the greatest interest for this study.

Petrography
Initial petrographic characterisation was performed under reflected light using a Nikon Eclipse LV100 Pol microscope, equipped with a Prior Proscan III motorised stage. Optical scans at a resolution of ~ 3 μm/px were recorded on all samples for further documentation using the same instrument. To identify unknown minerals and clarify paragenetic relationships, samples were then carbon coated and examined using an FEI Quanta 450F scanning electron microscope equipped with a Bruker Quantax EDX detector, housed at Adelaide Microscopy (University of Adelaide).

Electron-probe microanalysis
Sphalerite compositions were determined quantitatively using a Cameca SX-Five EPMA, equipped with five tunable wavelength-dispersive spectrometers, located at the University of Adelaide. The instrument ran the PeakSite v6.2 software for microscope operation, and the Probe for EPMA software (distributed by Probe Software Inc.) for all data acquisition and processing. Operating conditions were 20 kV/30 nA with a defocused beam diameter of 3 µm.
The full list of analysed elements along with primary and interference standards are given in Tables ESM1.2 to ESM1.4 in the Electronic supplementary material (ESM). Matrix corrections of Armstrong-Love/Scott φ(ρz) (Armstrong 1988) and Henke MACs were used for data reduction. Due to the complexity of off-peak interferences in sulfide minerals, all elements were acquired using a multipoint background fit, excepting Cu and Ni, which were acquired using a traditional 2-point linear fit.
Beam damage and element migration (e.g. for Cd, In) were monitored by using the Time Dependent Intensity (TDI) correction feature of Probe for EPMA (e.g. Donovan and Rowe 2005). The decay of X-ray counts over time was measured and modelled to return a t = 0 intercept, and from this a concentration could be calculated. Upon visual inspection, the X-ray counts did not appear to decay over time and thus no correction was applied.

Laser ablation-inductively coupled plasma-mass spectrometry
The trace-element contents of sphalerite, chalcopyrite and pyrite were analysed using an ESI NWR213 solid state laser, coupled to an Agilent 7900 ICP-MS at the University of Adelaide. The system used He as the carrier gas, which was mixed with Ar for plasma generation in the ICP-MS. Where possible, a minimum of 5 spot analyses were performed on each mineral generation in each sample. The following isotopes were monitored during each measurement (30 s background, 50 s ablation): 55 Mn, 57 Fe, 59 Co, 60 Ni, 63 Cu, 66 Zn, 69 Ga, 72 Ge, 73 Ge, 75 As, 77 Se, 107 Ag,111 Cd,113 In,115 In,118 Sn,121 Sb, 125 Te,202 Hg,205 Tl,208 Pb and 209 Bi. Ablation spot size was generally 50 μm. An ablation pulse frequency of 10 Hz and fluence of 3.5 J/mm 2 were used. In some samples, mineral grain sizes and intergrowth relationships required the use of smaller spot sizes down to 25 µm to avoid inclusions.
Data reduction was done using the Iolite software package (Woodhead et al. 2007;Paton et al. 2011). To convert measured count rates to concentrations, the 66 Zn signal was used as the internal standard for sphalerite, while 57 Fe was used for chalcopyrite and pyrite. Average Zn concentrations for each sphalerite generation in each sample were taken from EPMA measurements. Stoichiometric Fe concentrations were used for pyrite and chalcopyrite. MASS-1 (Wilson et al. 2002) was used as the external standard for all elements. In addition, two measurements on NIST SRM 610 (NIST 2012) were included with every standard block for quality control. A block of two to three standard measurements was inserted before and after every 20 to 30 sample measurements. Off-line corrections were made for the isobaric interferences of 113 Cd on 113 In and 115 Sn on 115 In in both the samples and the standard, using raw count rates and the natural abundance ratios of the relevant isotopes. Interferences of 56 Fe 16 O and 57 Fe 16 O (cf. Belissont et al. 2014) on 72 Ge and 73 Ge, respectively, were monitored by comparing the measured abundance ratio of 72 Ge to 73 Ge in the sample to the natural ratio expected for the two isotopes. If significant deviations (more than 50% relative) of 72 Ge/ 73 Ge from the natural value of 3.5 occurred, the measurement was designated as below detection limit, with the highest of the two reported concentration values as the detection limit. In general, interferences did not appear to be a problem for Ge measurements, as also noted by Belissont et al. (2014). Even in Fe-rich, Ge-poor sphalerites, spurious concentrations produced by the interference of Fe-O species on Ge never exceeded ~ 0.5 µg/g. The official reference value of 58 ppm Ge was used for the MASS-1 standard (cf. Belissont et al. 2014), and count rates on MASS-1 always showed the correct 72 Ge/ 73 Ge ratio within 10% relative of the expected natural value.
Since the Iolite software produces unrealistically low detection limits in cases where background counts are below the minimum count rate (50 cps in our case), median detection limits for the affected elements were estimated assuming a count detection limit of twice the minimum count ratei.e. 100 cps. This procedure was applied to the following elements: Co, Ni, Cu, Ga, Ge, Se, Ag, Cd, In, Sn, Sb, Tl and Pb.

Estimation of ore-forming conditions from sphalerite composition
Sphalerite formation temperatures were estimated using the GGIMFis geothermometer (Frenzel et al. 2016). This is based on an empirical relationship between sphalerite composition and formation temperature and is described by the following equation: with where ln(x) denotes the natural logarithm (base e), and c i is the concentration of trace element i in sphalerite, given in units of µg/g for Ga, Ge, Mn and In, and in wt.% for Fe. The relationship was calibrated using mean sphalerite compositions and microthermometric data for 51 hydrothermal base-metal sulfide deposits, formed between 100 and 400 °C (cf. Frenzel et al. 2016). Uncertainties on absolute temperatures estimated from (1) are generally on the order of ± 50 °C across the entire calibration range. The GGIMFis geothermometer is expected to work well for the estimation of average formation temperatures for individual deposits, as well as distinct mineralisation events within deposits (cf. Frenzel et al. 2016;Bauer et al. 2019). Its ability to capture smallerscale variations in formation temperatures-e.g., across an individual sphalerite grain-has not yet been tested.
In addition to formation temperatures, we also estimated the sulfur fugacities which prevailed during sphalerite formation. To do so, we first recalculated the measured Fe-contents in sphalerite to mol.% FeS using the Fe/Zn ratios determined by LA-ICP-MS. We then used these concentrations to estimate the FeS-activity in the sphalerite, FeS(sp) , from where FeS is the Fe-content in sphalerite (mol.% FeS), and FeS max is the temperature-and pressure-dependent maximum solubility of FeS in sphalerite when in contact with metallic Fe and troilite, described by: where T is temperature in Kelvin, and p is pressure in kbar ( Barton and Toulmin 1966;Balabin and Urusov 1995). A formation pressure of 0.2 ± 0.1 kbar (cf. Voudouris et al. 2008) and formation temperatures estimated from the GGIMFis geothermometer were used in Eq. (4). Finally, the resultant FeS-activities were combined with GGIMFis temperatures to calculate fS 2 values, including uncertainties, using: Equations (3) to (5) are based on a re-fitting of the data from Barton and Toulmin (1966), Scott and Barnes (1971), Scott and Kissin (1973), Balabin and Urusov (1995) and Lusk and Calder (2004) as described in detail in ESM 1.
This estimation procedure for fS 2 relies on the assumption that Fe contents in sphalerite are buffered by pyrite. Pyrite is present in most samples and is generally an abundant and widely distributed mineral at Plaka. It occurs throughout the complete paragenetic sequence (Voudouris et al. 2008). It is always present with Sp I in our samples, but not generally with Sp II (cf. "Results"). Nevertheless, its general abundance across the deposit means that the fluids involved in the precipitation of Sp II would have been in contact with pyrite on their way, and generally in close proximity to, the site of mineralisation, even if pyrite and sphalerite did not precipitate together. This means that pyrite-buffering is a reasonable assumption in the present case. The same working principle was adopted by Barton et al. (1977) for their classic study of the Creede district, Colorado.

Results
The following sections briefly summarise the key petrographic observations made on the investigated samples, as well as the corresponding results for sulfide mineral chemistry.

Petrography
Overall, the paragenetic sequence in all investigated samples is similar. Note, however, that not all parts of the sequence described below occur in each sample. For instance, sample PLA-7 only contains minerals of the late vein-related paragenesis, while sample 46 only contains minerals of the early carbonate-replacement paragenesis. Mineral names and/or abbreviations in combination with Roman numerals are used to indicate distinctive mineral generations occurring consistently across samples.
Quartz is generally the first hydrothermal mineral, occurring mostly as isolated crystals and clusters of euhedral grains. This is followed by an early pyrrhotite (Po I), now completely replaced by pyrite and/or marcasite ( Fig. 3B, C). An assemblage comprising Fe-rich sphalerite (Sp I) ± galena (Gn I) ± pyrrhotite (Po II) ± chalcopyrite (Ccp I) succeeds Po I, often in association with siderite ( Fig. 3B). Sphalerite I is generally characterised by abundant small inclusions of chalcopyrite ('disease', Barton and Bethke 1987, Fig. 3F), while Ccp I is characterised by small star-shaped inclusions of sphalerite (Fig. 3D). Pyrrhotite II is often partially replaced by Ccp I, as well as intimate intergrowths of sphalerite and pyrite (Fig. 3B). Euhedral to massive pyrite (Py I) concludes this early stage.
These observations are consistent with the descriptions of Voudouris et al. (2008) and confirm that the major mineralisation stages at Plaka are well covered by the investigated sample suite. Table 1 summarises the results of the sphalerite LA-ICP-MS analyses in terms of the geometric means for each sphalerite generation in each sample, accompanied by the lower and upper bounds of the corresponding 95% confidence intervals. These statistics were calculated by replacing all values below the detection limit (B.D.L.) with the respective detection limit (cf. van den Boogaart and Tolosana-Delgado, 2013) and omitting values affected by mineral inclusions.

Trace elements in sphalerite
Elements not shown in Table 1 were either mostly below or around detection limit (Ni, Se, Te, Tl, Bi) or were largely hosted by inclusions (As, Pb). A summary of EPMA results (Table ESM1.5) as well as files with the detailed documentation and analytical results for each EPMA/LA-ICP-MS measurement spot are included in the ESM.
It is apparent from Table 1 that there are systematic compositional differences between sphalerite I and sphalerite II. Where both generations occur in the same sample, concentrations of Fe, Mn, Co and In are generally higher in sphalerite I. Concentrations of Ga, Ge, Sn, Sb and Hg, on the other hand, show the opposite behaviour. Finally, Cu, Ag and Cd concentrations do not show clear systematic differences between the two sphalerite generations across samples. Another interesting observation is that high concentrations of Co and In appear to be restricted to sphalerite I in samples collected from the southern part of the deposit-i.e. in or close to Vein 80 (samples 38, FL8012, FL8015 and PLA-7; cf. Table ESM1.1).
Formation temperatures estimated from the GGIMFis geothermometer (Frenzel et al. 2016) are consistent with the chemical differences described above (Table 2). Namely, sphalerite I appears to have formed at a substantially higher temperature (276 ± 54 to 373 ± 61 °C on average) than sphalerite II (182 ± 50 to 212 ± 49 °C on average) in all samples. While relatively large uncertainties are attached to the estimates of absolute temperatures (Table 2), the differences between the two generations are well constrained, as shown by the substantial differences in PC1* values, which do not overlap at all between the generations, even across samples. This agrees well with previous observations from fluid inclusion studies which suggest the formation of sphalerite I is associated with early high-temperature fluids (~ 360 °C), while sphalerite II is associated with later low-temperature fluids (< 250 °C) (Voudouris et al. 2008).
Average sulfur fugacities estimated from GGIMFis temperatures and Fe contents in sphalerite (detailed procedures and assumptions in ESM1) also reveal some differences between the two sphalerite generations (Table 2). Specifically, sphalerite I seems to have formed at somewhat higher absolute sulfur fugacity (mean log 10 fS 2 values between − 11.9 ± 2.9 and − 8.1 ± 2.2) than sphalerite II (mean log 10 fS 2 values between − 16.3 ± 3.7 and − 14.2 ± 3.1). However, these differences are less pronounced than for the temperatures.  Table 1 Summary statistics of LA-ICP-MS data for sphalerite Mean, geometric mean; C.I. Min., lower bound of the 95% confidence interval of the mean; C.I. Max., upper bound of the 95% confidence interval of the mean. All statistics calculated from measurements unaffected by inclusions (smooth ablation traces); values below detection limit were replaced by the detection limit for the calculations (cf. van den Boogaart and Tolosana-Delgado, 2013) a Concentrations of Hg must be viewed as semi-quantitative due to large uncertainties regarding the Hg content of the MASS1 standard material. However, relative differences between concentrations should be reliable b Cu concentrations are likely too high for this sample due to abundant Ccp inclusions in this sphalerite generation, which could not be removed during data processing Sample Generation n Fe (wt.%) Mn (µg/g) Co (µg/g) Cu (µg/g) Ga (µg/g) Ge (µg/g) Ag (µg/g) Cd (µg/g) In (µg/g) Sn (µg/g) Sb (µg/g) Hg (µg/g) a

46
Sp I

Trace elements in chalcopyrite
A summary of the chalcopyrite data is presented in Table 3. Again, the data for elements that were largely around or below minimum detection limit (Mn, Ni, As, Sb, Te, Hg, Tl), or mostly hosted in inclusions (Pb), are not shown in this table. Complete results for each measurement spot are included in ESM2.
Similar to sphalerite, there are systematic chemical differences between the two generations of Table 2 GGIMFis temperatures and log 10 (fS 2 ) values calculated from sphalerite compositions (arithmetic means ± 2 σ) Note: Uncertainties provided here for PC1*, T GGIMFis , log 10 (FeS), and log 10 (fS 2 ) in the corresponding columns are internal uncertainties only, arising from variations in sphalerite composition within the generations within each sample. External uncertainties described by 2 σ external in separate columns contain all sources of uncertainty, including those arising from the calibration of the GGIMFis geothermometer (Frenzel et al. 2016) as well as the model describing log 10 (fS 2 ) as a function of sphalerite composition (cf. ESM 1). The uncertainty levels correspond to ± 2 standard errors in all cases. For log 10 (fS 2 ), most of the external uncertainty is in fact due to the uncertainty in absolute GGIMFis temperatures Sample

Trace elements in pyrite
The pyrite data are summarised in Table 4. Unlike sphalerite and chalcopyrite, only one generation of pyrite could be analysed in each sample. This was the early euhedral pyrite occurring slightly later than sphalerite I and chalcopyrite I in most samples. While it is apparent that pyrite compositions differ markedly from sample to sample, there are no clear systematic differences between samples taken from different parts of the deposit, except in terms of their Se contents. These are higher for the two samples taken in the south, in or around Vein 80. Furthermore, the only sample containing measurable Bi concentrations in pyrite is FL 8012, also from the south.

Discussion
Below, we briefly discuss the implications of our results for the evolution of the Plaka ore-forming system. We then focus on the LA-ICP-MS results for each mineral, how they relate to this evolution, and which implications this has for the overall behaviour of the different trace elements during ore formation.

Evolution of the Plaka ore-forming system
While the temperature and salinity evolution of the oreforming fluids at Plaka was previously constrained through microthermometric measurements by Voudouris et al. (2008), the additional determination of sulfur fugacities in this contribution (Table 2) enables a discussion in terms of the fS 2 -T paths of the fluids (Einaudi et al. 2003). Figure 4 shows the position of the different Plaka sphalerite generations in fS 2 -T space, relative to important mineral reaction lines. Despite the relatively high uncertainties in terms of their absolute formation temperatures, there are clear differences between the two sphalerite generations. While sphalerite I falls on the boundary between the lowsulfidation and intermediate-sulfidation fields, sphalerite II lies well within the intermediate-sulfidation field. Furthermore, the evolutionary trend defined by the two sphalerite generations follows the "rock buffer" of Einaudi et al. (2003) rather than the sulfur-gas buffer. This indicates that the ore fluids cooled in equilibrium with the magmatic and metamorphic country rocks, and that this was the main driver for the observed changes in ore-forming conditions (cf. Einaudi et al. 2003), in addition to mixing with meteoric and/or marine waters (cf. Voudouris et al 2008).

Sphalerite data
We showed earlier that mean sphalerite formation temperatures estimated from the GGIMFis geothermometer are consistent with the known evolution of the Plaka mineralising system, i.e., sphalerite I corresponds to the early hightemperature stage while sphalerite II corresponds to the later low-temperature stage. An interesting question is whether a Table 4 Summary statistics of LA-ICP-MS data for pyrite (in µg/g) Mean, geometric mean; C.I. Min., lower bound of the 95% confidence interval of the mean; C.I. Max., upper bound of the 95% confidence interval of the mean. All statistics calculated from measurements unaffected by inclusions (smooth ablation traces); values below detection limit were replaced by the detection limit for the calculations (cf. van  To explore this, Fig. 5 compares the distribution of GGIM-Fis temperatures obtained from individual LA-ICP-MS spot analyses to previously determined fluid-inclusion (FI) homogenisation temperatures in associated gangue minerals. It is evident that the GGIMFis temperatures reproduce the distribution of FI homogenisation temperatures, both in terms of the central tendency and spread of the high-and low-temperature populations. Given that sphalerite occurs throughout the entire paragenetic sequence of the deposit (cf. petrography section; Voudouris et al. 2008), and that LA-ICP-MS and microthermometric data are sampled at a similar scale (10 s of µm), this suggests that, indeed, sphalerite chemistry may record smallscale variations in sphalerite formation temperatures. This is an important result since the original analysis by Frenzel et al. (2016), as well as later work by Bauer et al. (2019), only showed that the GGIMFis geothermometer is applicable at the scale of entire ore deposits and/or temporally and physically very distinct mineralisation events within the same deposit.
We made further use of this result, in conjunction with the estimated sulfur fugacities for each measurement spot, to investigate the dominant controls on the observed variations in sphalerite composition. To get a first impression of the trends present in the dataset, we plotted trace-element concentrations for the individual measurement spots against temperature (Fig. 6) and a value describing variations in sulfur fugacity corrected for temperature-driven changes (Fig. 7). The correction of the temperature effect on sulfur fugacity before plotting was necessary since a substantial proportion of the total variance in log 10 (fS 2 ) values appear to be due to changes in temperature (Fig. 8A), probably in consequence of the cooling-driven evolution of the ore forming system (cf. previous section). This means it would not have been possible to distinguish between the specific effects of temperature and sulfur fugacity without an appropriate correction. The correction was done by fitting a model of the following form to the dataset: where D and E are constants. The residuals of this model correspond to the corrected dfS 2 values used in Fig. 7 and the further modelling. The fitted values for coefficients D and E are indicated in Fig. 8A. Figures 6 and 7 already indicate that some trace-element concentrations correlate strongly with temperature (e.g. Fe, Mn, Ga, Hg) and/or dfS 2 (e.g. Ge, In), while others show no such clear trends (e.g. Cd). However, it is also apparent that various complexities are present within the dataset. For instance, Co concentrations seem to correlate with temperature in samples from the southern part of the deposit (Fig. 6C, black, green and blue symbols), while samples from the northern part fall off this trend (red and yellow symbols). In fact, clear groupings are observed for individual samples for nearly all elements. These groupings illustrate the variable effects of differences between sample locations (cf. Dmitrijeva et al. 2018).
Separating the effects of temperature, dfS 2 , and sample location by standard correlation analysis is nearly impossible. Therefore, we used linear mixed-effects models (cf. Winter 2013) of the following form to quantify the different effects: where c i denotes the concentration of trace element i , T denotes the estimated GGIMFis temperature in Kelvin,dfS 2 denotes a corrected value for log 10 (fS 2 ) after de-correlation with 1∕T , and Sample is a categorical variable describing on which sample a measurement spot was obtained. The notation used in Eq. (7) corresponds to that used in the R (6) log 10 (fS 2 ) = D T + E (7) log 10 c i ∼ 1 T + dfS 2 + (1|Sample) Fig. 4 Sulfur fugacity-inverse temperature plot adapted from Einaudi et al. (2003) showing the location of the two sphalerite generations from Plaka relative to different mineral reaction lines (in black), the sulfur gas (S-gas) and rock buffers of Einaudi et al. (2003)  software suite (R Core Team 2017) to describe linear mixedeffects models. The same notation is employed in Dmitrijeva et al. (2018) and Godefroy-Rodriguez et al. (2020). Temperature and dfS 2 are treated as fixed effects, while Sample is treated as a random effect. Interaction terms between the different variables were not included in the modelling. The main reason for using linear mixed-effects models was, first, that the dataset is not balanced. That is, different numbers of observations are available from each sample and for each sphalerite generation (cf. Table 1). This can lead to distortions in the fitting of linear models unless it is properly accounted for (cf. Winter 2013;Dmitrijeva et al. 2018). The second reason was the necessity to account for random signals introduced into the dataset by variations between sample locations which are not related to variations in the fixed effects T and dfS 2 , and which would otherwise cause distortions in the results. For instance, the different Co concentrations between samples from the northern and southern parts of the deposit described above (cf. Figure 6C) are a manifestation of such effects. This was accounted for in the models by including Sample as a random effect (cf. discussions in Dmitrijeva et al. 2018;Godefroy-Rodriguez et al. 2020). Effectively, the Sample variable can be thought of as a proxy for the (temperature-and sulfur fugacity-corrected) location dependence of sphalerite compositions, giving a compound measure for the variance due to both larger (kilometre) and smaller (decimetre)-scale effects. Table 5 summarises the results of the modelling for each element, providing p-values for each of the fixed effects, as well as R 2 values for the fixed and random effects. Unfortunately, it is not mathematically possible to give separate estimates of R 2 for each of the fixed effects (cf. Godefroy- Rodriguez et al. 2020). A graphical representation of the overall goodness of fit for each element is provided in Fig. 9, where observed values are plotted versus predicted ones.
The results in Table 5 can be interpreted as follows: the p-values provide a measure of statistical significance: the smaller the p-value, the more probable the existence of a relationship with the respective explanatory variable. In the present case, p-values below 2 × 10 −3 were considered to indicate statistically significant relationships. The R 2 values, on the other hand, provide a measure of the proportion of the total variance explained by a given part of the model. Specifically, R 2 (1/T + dfS 2 ) corresponds to the proportion of the variance explained by temperature and fS 2 , while R 2 (Sample) corresponds to the proportion of the variance explained by the effects of sampling/locality. R 2 (Total) is the sum of A.D.L., above detection limit (i.e., measurable); B.D.L., below detection limit. Copper has been omitted from this and the following figures since it does not show clear relationships with any of the investigated explanatory variables (cf. Table 5) the two previous values and describes the goodness of fit for the overall model. Thus, the R 2 (Sample) value of 0.67 for Co means that differences between localities account for 67% of the observed variance in Co concentrations, while the R 2 (1/T + dfS 2 ) value of 0.27 means that temperature and sulfur fugacity account for 27%. The sum of these two values is 94% which corresponds to the proportion of the variance in Co concentrations explained by the overall model, i.e. R 2 (Total).
Considering first the R 2 values for the overall models, there are three broad groups of trace elements: (1) those whose variability is well explained by the models (R 2 > 0.80)-i.e. Mn, Fe, Co, Ga, Ge, In, Sb and Hg; (2) those whose variability is accounted for in large parts (R 2 > 0.60), while considerable uncertainties remain-i.e. Sn and Cd; (3) and, finally, those whose variability is only explained to some degree by the model (R 2 < 0.50), and for which the main cause(s) of the observed variability cannot currently be determined-i.e. Cu and Ag. Within these groups, great variability exists in terms of the dominant geological control(s), as indicated by the highly variable R 2 values for (1∕T + dfS 2 ) and Sample.
Several features of the results are worth noticing. First, it is of considerable interest that several elements (Mn, Co, Ge, In, Sb and maybe Hg) in addition to Fe show statistically significant relationships with dfS 2 . These effects have never been described before. Unfortunately, the degree to which dfS 2 influences trace-element concentrations can only be assessed qualitatively, since separate R 2 values cannot be calculated. In general, the relationships with dfS 2 appear to be weaker or less well constrained than those with T , as indicated by the larger p-values for dfS 2 .
Second, it is interesting that some temperature dependence is seen here for all trace elements, except Cu and Cd. This probably reflects the general importance of cooling in the evolution of the Plaka ore-forming system (cf. previous section).
Third, the results for Sample highlight the considerable effects of sample location on some trace-element concentrations, even after correction is made for the effects of temperature and sulfur fugacity. The most striking cases are Co, where this accounts for 67% of the observed variance, Cd, where it accounts for 65% and In, where it accounts for 61%. These effect sizes are similar in magnitude to those reported for hematite in the Middleback Ranges iron deposits, and pyrite in the Kalgoorlie gold district (Dmitrijeva et al. 2018;Godefroy-Rodriguez et al. 2020).
While the exact reasons for such location-dependent effects remain unclear, it is likely that they reflect local variations in the physico-chemical parameters of ore formation not accounted for in the modelling-e.g. pH, fO 2 and fluid composition. Potential explanations for such local differences include variations in host-rock composition, relative position with respect to the major hydrothermal feeder zones of the deposit, as well as compositional differences between locally introduced fluid pulses.
Finally, Table 6 shows the model coefficients for the dependence of different elements on T and dfS 2 . These provide an idea of the direction and degree in which changes in these parameters affect specific elements. Of particular interest here are the slopes with respect to dfS 2 for Mn, Co, Ge, In, Sb and Hg, since they provide direct indications for the relevant reaction mechanisms involved in their incorporation into the sphalerite (cf. ESM1 for the discussion on Fe). For instance, a slope of + 0.5 with dfS 2 indicates the consumption of around 0.5 S 2 molecules per atom of trace element incorporated into the sphalerite, while a slope of − 0.5 indicates the release of 0.5 S 2 molecules. These balances in turn indicate changes in the oxidation state of either the respective trace element or sulfur during the reaction.
Below, we give some tentative reactions which may explain the observed slopes for Mn, Co, Ge and Hg, and which would be compatible with the probable speciation(s) of these elements in the hydrothermal fluids, some host rocks and/or buffer minerals (Varekamp and Buseck 1984;Calvert and Pedersen 1996;Wood and Samson 2006;Liu et al. 2011), as well as their binding states in sphalerite (Grammatikopoulos et al. 2006;Cook et al. 2009;Bonnet et al. 2017) :   Fig. 7 Trace element concentrations in sphalerite as a function of corrected log-sulfur fugacity, dfS 2 , for individual LA-ICP-MS spot analyses: A Fe, B Mn, C In, D Ga, E Ge, F Cd, G Co, H Ag, I Hg, J Sb, and K Sn. Symbol sizes are larger than typical analytical uncertainties. A.D.L., above detection limit (i.e., measurable); B.D.L., below detection limit Note that the observed slope of Hg concentrations with dfS 2 is smaller than ~ 0.5, as indicated by reaction (11). This may be due to some of the Hg being present in solution as Hg 2+ (cf. Varekamp and Buseck 1984) which can be directly incorporated into sphalerite without any participation of S 2 in the reaction.
For indium, a coupled substitution reaction with Cu (cf. Johan 1988), buffered by chalcopyrite and pyrite may explain the positive sign of the slope but would predict a smaller absolute value of ~ 0.25: For Sb, however, reactions compatible with both the indicated slopes with dfS 2 and its speciation cannot presently be devised. Antimony is generally thought to occur with an oxidation state of + III in the hydrothermal fluids and sphalerite (Pokrovski et al. 2006;Wood and Samson 2006;Cook et al. 2012;Belissont et al. 2014). While native Sb or stibarsen could buffer its incorporation into sphalerite, the corresponding reaction should produce the opposite of the observed trend since it would consume S 2 . Therefore, the cause for the observed trend in Sb concentrations remains unclear.
For the other elements for which no statistically significant relationships with dfS 2 were found, it is probable that the oxidation states of these elements in the fluid and sphalerite are the same, that their incorporation is not associated with other elements whose oxidation state changes, or that their concentrations are not buffered by pyrite. Therefore, no consumption or release of S 2 is required for their incorporation into the sphalerite structure.
Last but not least, the slopes with 1∕T are also interesting since they indicate the direction in which temperature affects the different trace elements. Referring back to Table 6, concentrations of Mn, Fe, Co and In increase with temperature, while concentrations of Ga, Ge, Ag, Sn, Sb and Hg decrease. Concentrations of Cu and Cd are not affected. It is worth noting that the results for Cu, Cd, Ga, Ge, In, Mn and Fe are similar to those obtained for the global dataset of Frenzel et al. (2016). That is, the concentrations of Mn, Fe and In increase with temperature, those of Ga and Ge decrease, and those of Cu and Cd are not affected. Only Co and Ag differ from their respective global trends, showing some dependence on T in the present dataset, but none in the global data of Frenzel et al. (2016). This may be due to differences in the dominant geological controls at different scales of observation. For instance, differences in background signals (e.g. source rock composition) may be more important at the global scale and thus reduce the observable strength of the temperature signal when compared to observations at the local scale. Note that the strength of the local relationships of Co and Ag with T at Plaka is quite weak, with R 2 values of < 0.27 and < 0.23, respectively (cf. Table 5). This would support such an interpretation.

Chalcopyrite data
In general, chalcopyrite chemistry is not sufficiently well understood to enable a similarly detailed discussion as the sphalerite data. Particularly the dependence of trace-element concentrations on formation conditions has not yet been constrained, neither globally nor at the scale of individual ore deposits (George et al. , 2018. Therefore, we must restrict ourselves to (1) noting that all measured trace-element concentrations in chalcopyrite from Plaka fall well within the range reported for other magmatic-hydrothermal deposits (cf. Cook et al. 2011;George et al. 2016George et al. , 2018 and (2) comparing the observed trends in chalcopyrite compositions with those described for sphalerite above. Regarding this last point, the discussion below is separated into temporal and spatial trends. Temporal trends refer to the differences between mineral generations within individual samples, while spatial trends refer to systematic differences between samples, particularly those from the southern and northern parts of the deposit (cf. "Results").
In terms of temporal trends, it is interesting to note that some trace elements (Ge, In, maybe Ag) show the same trends in chalcopyrite as they do in sphalerite, while others show opposite trends (Co, Sn) (cf . Tables 1 and 3). A third group (Cd, Ga) shows dissimilar, but not opposite, behaviours.
While it is probable that the observed differences between the two chalcopyrite generations reflect a general control by fluid evolution, just as they do for sphalerite, it is not clear which exact factors cause the observed trends. Only one constraint can be derived from the present data: where the concentration of a specific trace element in chalcopyrite changes in the opposite direction compared to sphalerite during cooling, it is likely that a change in fluid composition-in particular, the relative activity of the trace element in question-is not the dominant cause for the observed behaviour, since this should have produced similar trends. The reason for this is that changes in the activity of a given trace-element in the fluid should always produce changes in the same direction for all minerals into which the element is incorporated (cf. McIntire 1963). For instance, if the activity of Co 2+ aq in the fluid increases, everything else staying the same, then its concentration in all co-existing minerals would be expected to increase. This response can be overruled, however, by simultaneous changes in partitioning behaviour, e.g. through temperature-or pressure-induced changes in the equilibrium partitioning constants, or due to kinetic effects (cf. McIntire 1963). Changes in the activity ratios of the major elements typically replaced by trace elements in the two minerals, i.e. Fe to Zn or Cu to Zn in the case of chalcopyrite and sphalerite, are another possibility to independently change the relative rates of incorporation for the different trace elements (cf. McIntire 1963;Frenzel et al. 2016 for detailed discussion).
Opposite behaviour between sphalerite and chalcopyrite at Plaka is only observed for Sn and Co. Therefore, it is likely that changes in the activities of Sn and Co in the fluids between the earlier and later mineralisation phases are not the main driver for the observed temporal behaviour of these two elements. Instead, this is likely to be driven by other factors such as temperature changes.
For the other trace elements, no further constraints can be derived from the observed trends. Concentration changes in the same direction can be due to any number of factors, from changes in the activity of the relevant trace element itself, to changes in temperature, Fe, Cu and Zn activities or kinetic factors influencing trace-element incorporation. Clearly, there is a need for further work to improve the general understanding of chalcopyrite chemistry to better constrain the causes of the observed temporal trends.
In terms of spatial trends, it is interesting that some trace elements behave similarly in coexisting sphalerite and chalcopyrite. Particularly, Co and In show systematically higher concentrations in samples from the southern part of the deposit (cf. "Results"). The most probable explanation for these trends are differences in fluid composition between the locations. That is, fluids in the southern part of the deposit may have had higher activities of Co and In. In turn, these higher activities may reflect either differences in background signals (e.g. host rock compositions) or fluid pathways (e.g. proximity to fluid source, interaction with different rocks during ascent) between the locations. Finally, it is interesting to note in this context, that the partitioning of In between coexisting chalcopyrite and sphalerite at Plaka is similar to that reported from other deposits Carvalho et al. 2018;Frenzel et al. 2019). That is, sphalerite, with one exception, contains between 1 and 10 times higher In concentrations than coexisting chalcopyrite (cf. Tables 1 and 3).

Pyrite data
As for chalcopyrite, the current understanding of the physicochemical controls on pyrite composition in hydrothermal systems remains limited. Despite the relatively large amounts of analytical data collected over the past 30 years (e.g. Bajwah et al. 1987;Reich et al. 2005;Large et al. 2009;Genna and Gaboury 2015;Steadman et al. 2021), no meta-analysis has ever been conducted to systematically constrain the differences between different deposit types or the relationships with important physical parameters such as formation temperature. Furthermore, there is no comprehensive quantitative understanding of the thermodynamic controls for any of the common trace elements in pyrite such as exists, for instance, for Fe in sphalerite. Experimental data are only available for As and Au incorporation at low temperatures (Kusebauch et al. 2018(Kusebauch et al. , 2019. Theoretical predictions exist for Se at constant fluid composition (Huston et al. 1995). Broad correlations between formation temperatures and some trace elements (As, Au, Se) in natural pyrites have also been reported, but only from a limited number of Cu-Au deposits (Deditius et al. 2014;Keith et al. 2018). At present, it is not clear whether these relationships apply to pyrite in Pb-Zn deposits.
Given these limitations, we must restrict ourselves to some general comments on the Plaka data. First, we note that the observed As and Se concentrations in pyrite I would be compatible with formation temperatures above 300 °C, if indeed the relationships reported by Deditius et al. (2014) and Keith et al. (2018) apply here. This would be consistent with the other temperature data (GGIMFis, FI thermometry) available for the early high-temperature phase with which pyrite I is associated.
Second, the average Co/Ni ratio of pyrite I varies between 1 and 10 (Table 4). This fits the magmatichydrothermal origin of the ores (Bajwah et al. 1987;Gregory et al. 2015).
Finally, it is worth noting that pyrite composition also shows some systematic differences between the southern and northern parts of the deposit, some of which coincide with the trends already described for sphalerite and chalcopyrite. Specifically, higher Co concentrations in both sphalerite and chalcopyrite from the southern part of the deposit are also reflected in pyrite. Higher Se concentrations in pyrite from the south mirror the trend observed in chalcopyrite. Again, these features probably reflect systematic differences in the activities of the respective elements in the ore-forming fluids, with similar causes as discussed for chalcopyrite above.

Conclusions
Overall, this article has highlighted the great scientific potential which studies of sulfide mineral chemistry have in constraining both the evolution of hydrothermal oreforming systems, as well as the behaviour of different trace elements within them. Not only was it possible to reconstruct the T-fS 2 evolution of the Plaka system using sphalerite chemistry, with results that are consistent with previous fluid-inclusion studies, but it was also possible to gain important insights into the geological controls on the temporal and spatial distribution of several trace elements within the ore-forming system. Some of the results, such as the fS 2 dependence of Mn, Co, Ge and Hg concentrations in sphalerite, even provided new insights into the potential incorporation mechanisms of these elements into sulfide minerals.
Some of the insights reported in this paper, particularly with respect to the spatial distribution of trace elements within the deposit, could only be achieved by combining data from several minerals and comparing them in terms of the observed spatial and temporal trends. This demonstrates that there are some clear advantages in multi-mineral studies compared to the single-mineral ones which have mostly been the standard in the field so far (e.g. Cook et al. 2009;George et al. 2015George et al. , 2018Keith et al. 2018). Particularly, the combination of data on well understood minerals, such as sphalerite, with data on less well-understood minerals, such as chalcopyrite and pyrite, may prove a powerful approach to yield important new insights into ore-forming processes in the future.
Funding Open Access funding enabled and organized by Projekt DEAL. M. Frenzel's position was funded by the German Academic Exchange Service (DAAD) for the duration of his stay in Adelaide in 2017/18.

Data availability
All analytical data are included with this article in electronic supplements. Sample material used for this article is stored by M. Frenzel and P. Voudouris and can be made available on request.
Code availability No programs were written for this contribution. Fig. 9 Comparison between fitted and observed trace-element concentrations in sphalerite, according to the model described in the main text (Eq. (7); cf. Tables 5 and 6): A Mn, B Co, C Ga, D Ge, E Ag, F Cd, G In, H Sn, I Sb, and J Hg. Indicated R 2 values correspond to R 2 (total) in Table 5. A.D.L., above detection limit (i.e., measurable); B.D.L., below detection limit. Note that displaced stacks of values for individual samples (e.g. F) essentially reflect the effect of Sample in the models. Also note that the true relationship between Ge concentrations and the explanatory variables considered here is probably stronger, since the relatively large proportion of B.D.L. values tends to obscure part of the relationship. Iron has been omitted from this figure since its behaviour, by definition, is perfectly explained by temperature and sulfur fugacity ◂ Table 6 Fitted model coefficients for T and dfS 2 Note: A and B correspond to the best fits of slopes of the data with the respective fixed effect ± 2 standard errors. Negative slopes with 1/T indicate a positive relationship with T, and vice versa