Quasiparticle Interference and Symmetry of Superconducting Order Parameter in Strongly Electron-Doped Iron-based Superconductors

Motivated by recent experimental reports of significant spin-orbit coupling (SOC) and a sign-changing order-parameter in the Li$_{1-x}$Fe$_x$(OHFe)$_{1-y}$Zn$_y$Se superconductor with only electron pockets present, we study the possible Cooper-pairing symmetries and their quasiparticle interference (QPI) signatures. We find that each of the resulting states - $s$-wave, $d$-wave and helical $p$-wave - can have a fully gapped density of states (DOS) consistent with angle-resolved photoemission spectroscopy (ARPES) experiments and, due to spin-orbit coupling, are a mixture of spin singlet and triplet components leading to intra- and inter-band features in the QPI signal. Analyzing predicted QPI patterns we find that only the spin-triplet dominated even parity $A_{1g}$ (s-wave) and $B_{2g}$ (d-wave) pairing states are consistent with the experimental data. Additionally, we show that these states can indeed be realized in a microscopic model with atomic-like interactions and study their possible signatures in spin-resolved STM experiments.


I. INTRODUCTION
In iron-based superconductors, it has been widely believed that superconductivity is driven by repulsive interactions, enhanced by the presence of the spin fluctuations associated with the parent antiferromagnetic state. In this scenario, these fluctuations drive a sign reversal (s ± state) 1-5 between order parameters (OP) on the electron and hole Fermi surface pockets at the M -and Γ-point, respectively. The discovery of superconductivity in intercalated or monolayer FeSe at a critical temperature of the order above 40K revived interest in Fe-based superconductivity, but raised further questions on the origin of superconductivity in these compounds [6][7][8][9][10] , because, unlike bulk FeSe, ARPES experiments show that many of these FeSe-derived systems appear to be missing the hole pockets at the Γ-point required in the the conventional scenario.
Initial model calculations based on the multiorbital spin-fluctuation framework for systems manifesting only electron pockets at the M -point predicted d-wave symmetry state in this case 11,12 , driven by the spin fluctuations connecting the electron pockets that remain when the hole pockets are removed. In the proper 2-Fe unit cell, such a state must have gap nodes on the Fermi surface 13 . This is because the electron pockets located near (π, 0) or (0, π) points of the Brillouin Zone (BZ) in the 1-Fe unit cell fold onto (π, π) point of the folded BZ as the crystallographic symmetry lowers due to the Se positions. This may lead to hybridization between the electron pockets 14 , which then forces the d x 2 −y 2 -state to acquire gap nodes, although in principle the nodal area may be very small, proportional to the hybridization ("quasinodes"). On the other hand, ARPES experiments in most of the electron-intercalated materials indicated a nodeless superconducting (sc) state 15,16 . Several proposals for the gap structure have been put forward, including a conventional s ++ -wave scenario based on the electron-phonon interaction and orbital fluctuations 17 , as well as the "bonding-antibonding" scenario 13,18 in which the order parameter on the inner electron pocket (mostly d xz /d yz character) has one sign, and on the outer electron pocket (mostly d xy character) the other 14 . Furthermore, it has been argued that the hybridization of the electron pockets is mainly due to SOC [19][20][21] , which within a 3D spin fluctuation framework may stabilize the bonding-antibonding s ± state against d-wave 19 and mixes a spin-triplet component into the even parity s +− -wave state 20,21 . Overall, one can see that the sign structure of the superconducting order parameter is intimately related to the pairing mechanism. Therefore, experiments allowing to determine it could be of great potential importance.
One rapidly developing technique to determine the phase structure of the order parameter makes use of QPI as measured by Fourier transform scanning tunneling microscopy (FT-STM). This probe measures the wavelengths of Friedel oscillations caused by impurities present in a metallic or superconducting system, which in turn contains information on the electronic structure of the pure system. A subset of scattering wave vectors q corresponding to peaks in the FT-STM can be enhanced or not according to the type of disorder and the phase structure of the superconducting gap 22,23 . Recently it was proposed by Hirschfeld, Altenfeld, Eremin and Mazin (HAEM) 24 that the sign structure of the order parameter in a multiband system can be extracted from the Fourier transform of the local density of states QPI pattern near an impurity in the superconducting state. The antisymmetrized QPI intensity integrated over the wavevectors corresponding to scattering between two bands was shown to have a dependence on frequency very different for sign-changing and sign preserving scenarios leading to a strong, single-sign enhancement of the integrated response in the former case. This qualitative result was also confirmed by extensive numerical simulations with finite disorder 25 . Recently, a complementary phase sensitive technique to detect sign-changing gaps in the presence of strong impurity bound states was proposed 26 .
Using QPI analysis, the authors of Ref. 27 were able to identify a sign changing order parameter in FeSe. Most importantly for our purposes, similar conclusions were recently drawn for the strongly electron doped iron-based superconductor lithium hydroxide intercalated FeSe 28 . In other words, the order parameter in Li 1−x Fe x (OHFe) 1−y Zn y Se, alternates sign, either between the Fermi surface sheets, or within individual sheets. However, distinguishing between these alternatives was beyond the resolution of the experiment. In any case, the situation is somewhat more complicated than anticipated in Ref. 28, since the effect of spin-orbit interaction on pairing needs to be taken into account as well. Moreover, recent observation of Majorana zero modes in (Li 0.84 Fe 0.16 )OHFeSe 29 suggests possible broader implications of the spin-orbit coupling for the Cooper-pairing in electron doped intercalated iron-based superconductors. Note that in contrast to Ref. 28, no Zn substitution was used in Ref. 29. The amount of Zn, however, is relatively small (less than 2 percent). This amount does not affect the superconducting transition temperature or the electronic structure in a significant way and is done only for the purpose of enhancing the QPI signal in the scanning tunneling microscopy.
In this manuscript we study the possible Cooper-pairing symmetries and their QPI signatures for strongly electrondoped Fe-based superconductors using the effective three-orbital model of Refs. 20 and 21 with spin-orbit coupling and proper consideration of all lattice symmetries of the FeSe space group. We find that each of the resulting states -A 1g -wave, B 2g -wave and helical E u -wave -can have a fully gapped DOS consistent with ARPES experiments and, due to spin-orbit coupling, are a mixture of spin singlet and triplet components leading to intra-and inter-band features in the QPI signal. Analyzing predicted QPI patterns we find that A 1g -wave pairing state, with the two dominant peaks in the DOS roughly corresponding to the gap energies on each pocket, and B 2g -wave pairing state both with a significant even parity spin triplet component are consistent with the experimental data. Moreover, we show that pairing states with dominant spin triplet component can be identified using spin-resolved STM. 1. (a) Single layer of the iron based superconductors lattice structure. Red and green dots are iron and pnictogen or Se atoms, respectively. One pnictogen sublattice is puckered above the iron layer (filled green dots) one is puckered bellow (empty green dots) which divides the iron atoms into sublattices A and B. One-iron unit cells for sublattices A and B are denoted by dashed squares. The two-iron unit cell, taking the puckering into account, is shown by the solid square. The vector τ0 connects sublattices A and B. We wish to describe the low energy states near the M,-points of the Brillouin zone using the orbitally projected band model of Ref. 20 and 21 for the two-iron unit cell. Near the Fermi level, only the xz, yz and xy orbitals contribute significantly; hence, the full 10-orbital tight binding model is projected onto the subspace of these three orbitals. The effective low energy Hamiltonian near the M-point that takes into account all the lattice symmetries of the FeSe space group as well as time reversal symmetry is defined as

II. MODEL
where the four component spinor Ψ † M,σ (k) = Ψ † X,σ (k), Ψ † Y,σ (k) describes the states at the M-point for each spin projection σ. The doublets Ψ X,σ (k) and Ψ Y,σ (k) are defined as Moreover, we have and where the Pauli matrices {σ x , σ y , σ z } and {τ 1 , τ 2 , τ 3 } act on spin and orbital space, respectively. The λ z , p z1 and p z2 terms in eq.(3) describe the k-dependent intra-band SOC which does not couple the two Fermi pockets but lifts the out-of plane spin degeneracy. The inter-band SOC term which hybridizes X-and Y-pocket is given by In order to describe intercalated FeSe we use the Luttinger invariants, Tab.(I), which were evaluated in Ref. 21 based on the available ARPES data.  The value of the intraband spin-orbit coupling, λ z is in agreement with those found in ab-initio calculations and ARPES experiments 30,31 . Furthermore, the value of the interband spin-orbit coupling, λ SOC between the electron pockets separated by the large momentum yielding their hybridization and splitting on the Fermi surface is found to be smaller and is taken to be ∼ 5 meV 31 .
Diagonalizing Eq.(1) yields four bands: two regular ones that form the inner and outer electron pocket, see Fig.1(b), and two incipient bands that do not cross the Fermi level. The effect of inter-band SOC on the band structure is visualized in Fig.1(c) and Fig.1(d).

III. MEAN FIELD PHASE DIAGRAM
Although phenomenologically the classification of superconducting orders for two electron pockets was considered previously 20,21 we analyze here its microscopic formulation via mean-field treatment of the atomic on-site interactions given by the Hubbard and Hund's couplings U , U ′ , J and J ′ which enter the Hubbard-Kanamori Hamiltonian as Here α ∈ {A, B}, {σ, σ ′ } ∈ {↑, ↓} and {µ, ν} ∈ {yz, xz, xy} label lattice sites, spins and orbitals, respectively. d α † µ,σ (r) and d α µ,σ (r) are the second quantized operators creating and annihilating particles on sub-lattice A and B, see Fig.1(a). Using the results presented in Ref. 20 and assuming sharply localized Wannier functions of the xz, yz and xy orbitals, we can relate, up to a constant, the d α µ,σ (k) operators acting in the one iron unit cell to the components of the doublets Ψ X,σ (k) and Ψ Y,σ (k) in the two iron unit cell via where we absorb the constant prefactors into the Hubbard and Hund terms. As the xy orbital contributes to both Xand Y-point eq.(10) leads to "Umklapp" terms at the M-point. We assume that even parity solutions are still the leading pairing instabilities (for odd parity solutions see Appendix B) and use Eqs. (8)(9)(10) to project eq.(7) onto the low energy model decoupled into the spin singlet A 1g s-wave and B 2g d-wave symmetry states which at the M-point in presence of inter-band SOC couple to the E g even parity spin triplet state. When defining the two doublets Ψ T 1σ (k) = (c yzσ (k), c xzσ (k)) and Ψ T 3σ (k) = (c xy X σ (k), c xy Y σ (k)) the pairing terms read The 2x2 block in eq.(11) corresponds to A 1g s-wave pairing. We write J ′ 13 = αJ ′ 11 and find two eigenvalues E A1g = 1 4 J ′ 11 2 + 3U ± (J ′ 11 ) 2 + α16J ′ 11 − 2J ′ 11 U + U 2 which correspond to ordinary "plus-plus" (s ++ ) s-wave and sign-changing "plus-minus" (s ± ) s-wave pairing, respectively. While the former channel is purely repulsive without spin-orbit coupling the latter becomes attractive once J ′ 11 > (U + U √ 1 + 8α)/4α. The paring term that leads to B 2g d-wave is E B2g = 1 2 (U − J ′ 11 ) and can be directly read off. It is attractive once J ′ 11 > U and competes with sign-changing s-wave. Since we assume sharply located Wannier functions which yields eqs.(8-10) and on-site interactions only we find that within our simple mean field approximation the xy orbitals do not contribute to d-wave pairing as no "pair-hopping" term J ′ 33 mediates between xy X and xy Y . This changes once the higher-order diagrams (spin fluctuations) are taken into account.
The E g even parity spin triplet corresponds to pairing between the first (second) component of Ψ Xσ (k) (Ψ Y σ (k)) and the second (first) component of Ψ Y σ (k) (Ψ Xσ (k)) and is thus inter-band and attractive once E Eg = U ′ 13 −J 13 < 0.
We perform a mean-field decoupling of Eq.(11) into A 1g and B 2g spin singlet channels with the pairing terms given by In terms of Ψ M,σ (k) a triplet term can be written as Ψ T M,σ (−k)Miσ y σΨ M,σ (k), whereM and iσ y σ represent orbital and spin part, respectively. Since iσ y σ is symmetric an even (odd) parity triplet requiresM to be antisymmetric (symmetric). iσ y σ can be divided into an in plane iσ y (σ x , σ y ) and out of plane iσ y σ z component which transform as as the two dimensional E g and one dimensional A 2g irreducible representation, respectively. In presence of SOC, orbital and spin degrees of freedom transform together under operations of the space group. We focus on the E g even parity spin triplet that together with the E g in-plane spin component decomposes into a direct sum of one dimensional representations as Using the two anti symmetric components of E g and τ ± = (τ 1 ± iτ 2 )/2 one finds two even parity spin triplets that transform according to A 1g and B 2g and hence, couple to the singlet channel.
We refer the reader to Refs. 21,32 for further details.
In terms of the spinor Ψ T with and the pairing terms∆ where the gaps in orbital space are given by the equations (A1) and (A2). We self-consistently compute∆ A1g and ∆ B2g as a function of temperature and inter-band SOC for the two cases E A1g < E B2g and E A1g > E B2g and present In the case when solely E A1g (E B2g ) is attractive, ∆ A t (∆ B t ) is induced by SOC and scales with λ SOC . In this case |∆ | and we call the state a spin singlet dominated A 1g (B 2g ). If E Eg is attractive, however, ∆ t can develop independently of ∆ 1,3 and allows for states where |∆ . Equally whether the system is in the singlet or triplet dominated regime, the form of Eqs. (21) and (22) stays the same.
Depending on the ratio α = J ′ 11 /J ′ 13 the leading pairing symmetry for low values of λ SOC is either B 2g d-wave or A 1g s-wave with accidental nodes and a sign change between inner and outer Fermi pocket. In Fig.2(a) we show a SOC mediated transition from B 2g d-wave to A 1g s-wave. For small values of λ SOC we find a nodeless d-wave state where the nodes are lifted and the gap is opened by a combination of spin singlet and triplet components with the singlet gap being dominant |∆ B 1 | > |∆ B t |. With increasing λ SOC the paring symmetry changes from nodeless to nodal d-wave and finally to nodeless s-wave with a dominant inter-band triplet component Fig.2(b) we choose parameters such that the initial state at small λ SOC is A 1g s-wave and has a dominant singlet gap Increasing λ SOC lifts the nodes at the Fermi level and drives the system in an A 1g s-wave symmetry with a dominant triplet gap |∆ A t | > |∆ A 1 |, |∆ A 3 |. We find that for large inter-band SOC the triplet dominated s-wave state wins over triplet dominated d-wave. The presence of an intra-band SOC term λ z does not change the phase diagram qualitatively. However, in the region where λ SOC is strong superconductivity in the d-wave channel is suppressed if λ z > 0 which further stabilizes the A 1g -wave solution. Moreover, a large triplet gap, in addition to SOC, lifts accidental nodes leading to a nodeless s ± pairing symmetry. Figure2(a) is consistent with the phase diagram of phenomenological model in Ref. 14 with several important differences. In particular, in our model the singlet dominated d-wave state competes with singlet dominated bonding anti-bonding s ± which arises from pair hopping between xz(yz) and xy orbitals. Furthermore, with increasing spinorbit coupling strength the A 1g state contains dominant interband spin triplet component which was absent in the simplified analysis of Ref. 14 as its pairing channel turned out to be strongly repulsive in the two band model. Since SOC couples even parity inter-band triplet to intra-band even parity spin singlet pairing, an attraction in the former induces a gap in the latter channel and vice versa. Consequently a gap at the FS opens even though the E g triplet is inter-band. A similar paring state with attraction in the triplet channel was recently proposed for highly doped systems with only hole pockets 33 .

A. Local density of states in presence of impurities
In order to investigate the QPI signatures of the possible pairing states we need to calculate the local density of states (LDOS) for a multi-orbital system in the superconducting state. We further investigate corrections to the LDOS which arise from scattering at single charge impurities to make a statement whether there is a sign-change between inner and outer electron pocket using HAEM's method. Moreover, we consider scattering from a single magnetic impurity to reveal information about the spin structure in the system.
In order to describe superconductivity and spin-resolved STM we introduce the 16 component Balian-Werthammer spinor Within this basis the Hamiltonian in the superconducting state is given by where the additional factor of 1 2 accounts for double counting and From Eq.(24) we find the superconducting Green's functionĜ 0 ] −1 and the local density of states (LDOS) given by where we assumed sharply localized Wannier functions of the xz, yz and xy orbitals. We now introduce a single on-site non-magnetic potential scatterer to the system which is located on either sublattice A or B. In the one-iron unit cell it can be described as We assume that the major contribution to scattering is of intra-orbital nature and project eq.(27) onto the states in the two iron unit cell. Due to the fact that d A/B xy in the projected model is a linear combination of c xy X and c xy Y the intra-orbital term induces scattering between xy X and xy Y thus also contributes to inter-band scattering. In the two iron unit cell the impurity potential is given by Following the HAEM 24 approach, we compute the antisymmetrized correction to the LDOS due to impurity scattering with δĜ q (ω) = kĜ 0 k (ω)ÛĜ 0 k+q (ω) being the convolution of the bare Green's functions dressed by a Nambu scattering matrixÛ = τ 3 ⊗ σ 0 ⊗Ṽ in Born approximation, i.e. with V α,β ≪ E F , where E F is the Fermi energy.
The HAEM method states that the momentum integrated and antisymmetrized LDOS, δρ − (ω), qualitatively differs in case of a sign-changing OP from that of a sign-preserving one yielding a strong enhancement of the integrated response for the sign-changing but not the sign-preserving one. This method has been recently successfully applied to confirm the sign-changing nature of the order parameter in FeSe 27 , where the superconducting gaps are also extremely anisotropic. More recently, the method has also been applied to (Li 1−x Fe x )OHFe 1−y Zn y Se 28 where δρ − (ω) shows a strong signal and no sign-change between 8 meV and 14 meV suggesting a sign-changing s ± pairing symmetry. Since our previous analysis shows that most of the conclusions regarding the phase structure of the superconducting gap obtained within Born limit are robust and remain valid also well beyond this limit, we restrict our analysis to weak potential scatterers.

B. Local density of states and phase sensitive correction to QPI
We would like to discuss singlet and triplet dominated A 1g -and B 2g -wave with respect to their consistency with QPI experiments in (Li 1−x Fe x )OHFe 1−y Zn y Se. Recently in Ref. (21) possible pairing states for intercalated FeSe have been discussed. Following Ref. [21] the pairing has to obey three criteria in order to be consistent with experiments: i) fully gaped LDOS, ii) the quasiparticle energy extrema are at or close to k F of the normal state ("back bending") and iii) two peak features in the DOS with a peak at about 8 meV and 14 meV, respectively. In the context of monolayer and intercalated FeSe it has been suggested that the second peak in the DOS arises due to pure interband gap 21 as in the former the gap size seen by ARPES 16 (13.7 meV) significantly deviates from the peak energy (20.1 meV) seen by STM 8 . This, however, does not have to be the case with the electron-intercalated materials where a gap of 13 ± 2 meV around the Fermi energy was reported 10 . Three pairing states were proposed to be consistent with experimental data. These are A 1g s-wave and B 2g d-wave and E u ⊗ U (1) helical p-wave. The A 1g and B 2g state were assumed to have a dominant intra-band singlet gap in order to ensure a full gap and back bending. An inter-band SOC then mixes in inter-band triplet E g pairing. For the odd parity intra-band triplet E u ⊗ U (1) p-wave state its the odd parity inter-band singlet A 2u that is coupled via SOC.
In the following we investigate spin singlet and spin triplet dominated A 1g -and B 2g -wave for their consistency with the citeria i)-iii) and the QPI data. For that we use Eq.(24) and calculate the LDOS (ρ(ω)) and the superconducting band dispersion to ensure i)-iii) are fulfilled. In addition, for each pairing state, we calculate the antisymmetrized correction to the LDOS (δρ − (ω)) to check whether the order parameter changes sign between electron pockets. Since the gaps in the intercalated FeSe are found to be isotropic 10 we further present the angular dependence of the gap projected on the inner and outer electron pocket, respectively. In Appendix (B) we also briefly discuss the QPI data for the odd parity E u -wave state and show that at least within Born scattering the results are not compatible with experiment. We find that A 1g and B 2g states with a dominant spin triplet gap can show back bending and, in contrast to their singlet dominated states and the odd parity E u -state, give an appropriate description of the QPI data found in experiment.
One should also further note that for the HAEM method the crucial role is played by the gaps present on the Fermi surface pockets. In the present case there is also additional interband gap. We show in the Appendix C that its phase structure with respect to the gaps present on the Fermi level cannot be elucidated within phase-sensitive QPI analysis. Furthermore, we show that the sign-changing and sign-preserving gaps on the Fermi surface still determine the characteristic features of δρ − (ω).

C. A1g-symmetry state
We start by examining spin singlet and triplet dominated s-wave pairing state. We use the singlet gaps ∆ 1 and ∆ 3 and the triplet gap ∆ t to fit band-structure and LDOS. To start with singlet dominated s-wave we utilize the fitting parameters presented in Ref. 21. As shown in Fig.3(a), ρ(ω) is indeed fully gaped and has a peak at 8 meV and 14 meV. The first peak between 6 − 10 meV mostly comes from the intra-band gaps at the electron pockets as can be seen by comparison with the gap projections on the inner, ∆ in and outer, ∆ out , electron pocket, Fig.3(c). At θ F = π/4 both gaps are equal in magnitude leading to a peak at ω ∼ 8 meV in ρ(ω). The width of the peak is limited by 6 meV and 10 meV which correspond to the minimum of ∆ in and the maximum of ∆ out , respectively, so that the occupied states below and above 8 meV are due to anisotropy of the gaps. The second peak at 14 meV is indeed due to inter-band pairing. It vanishes for λ SOC → 0 and ∆ t → 0. We find that both gaps are highly anisotropic and most important positive on both Fermi surface pockets, yielding an s ++ -pairing symmetry. In Fig.3(b) δρ − (ω) exhibits a sign change in the region 6 meV < ω < 10 meV which can be interpreted as a sign preserving order parameter and thus is not compatible with the QPI data 28 . Note that this state does not appear in our phase diagram, shown in Fig.2, where the A 1g state is sign-changing. The feature at ω ∼ 14 meV in δρ − (ω) for this state appears due to the small inter-band contribution and does not carry phase information as we explain for a simple model in Appendix C. In Fig.3(d) the superconducting band dispersion in ΓM -direction (θ F = π/4) is shown where the minima of the lower superconducting band are located above the former Fermi level (back-bending).
It turns out that the fitting parameters for spin singlet dominated s-wave, proposed in Ref. 21 would give rise to s ++ superconductivity. A possible s ± pairing symmetry with a dominant spin singlet gap, on the contrary, would require sign∆ A 1 = sign∆ A 3 . This state in our anaylsis exhibits highly anisotropic gaps and possible accidental nodes (see for example Fig.2c(1)). This makes it difficult to fit the U-shaped two peaked LDOS in the LiOH-intercalated FeSe data without invoking a spin triplet component, induced by spin-orbit coupling that lifts the nodes.

Spin-triplet dominated A1g-wave state
In Fig.4(a) and Fig.4(d) the LDOS and the superconducting band dispersion are shown. The LDOS shows a full gap and a peak at 8 meV and 13.7 meV, respectively, while the lower superconducting band shows back bending. The band projected gaps ∆ in and ∆ out we plot in Fig.4(c). Possible accidental nodes are lifted by a large ∆ t leading to a nodeless s ± pairing symmetry with two almost isotropic gaps. In Fig.4(b) we present δρ − (ω) which in contrast to the singlet dominated case does not change sign between 8 meV and 14 meV. If we compare Fig.4(a) and Fig.4(c) we find that the peak positions of the two peaks roughly agree with the the magnitudes of ∆ in and ∆ out and thus δρ − (ω) and the behavior of δρ − (ω) agrees with the experimental one. The first peak is sharper than the second since ∆ in is very isotropic. The second peak is broader as a consequence of an anisotropic ∆ out and the contribution from inter-band pairing. In contrast to the singlet dominated case δρ − (ω) exhibits well pronounced negative peaks at ω 1 ≈ ∆ out and ω 2 ≈ ∆ in without a sign change between them indicating an OP that changes sign between the electron pockets consistent with the experimental QPI data and ARPES data 10,28 . Also note that for the spin triplet dominated A 1g state the direct gap feature caused by ∆ out is very close to the interband feature in the DOS which agrees with the findings in Ref.10 and 28. The lower orange band has its minima centered above the former Fermi-level "back-bending".

D. B2g-symmetry state
We now discuss the singlet and triplet dominated B 2g d-wave pairing state. We use the singlet gaps ∆ 1 and ∆ 3 and the triplet gap ∆ t to fit band structure and LDOS. Note that even though our simple mean field approach leads to ∆ 3 = 0 for the d-wave case, as it involves decoupling of the repulsive intraorbital Hubbard interaction, a nonzero ∆ 3 is necessary to achieve qualitative agreement when fitting the experimental data. In theoretical calculations, a nonzero ∆ 3 appears due to inclusion of the spin fluctuation diagrams in the Cooper-pairing channel yielding momentum dependent interaction, see Ref. 11.

Spin-singlet dominated B2g-state
For the spin singlet dominated B 2g -wave state appropriate fitting parameters were found in Ref. (21). As shown in Fig.5(a) and Fig.5(d) the LDOS is fully gaped and has a peak at 8 meV and 14 meV, moreover, the lower superconducting band shows back bending. ∆ in and ∆ out are plotted in Fig.5(c). Both gaps are almost equal in magnitude and nodeless due to smallness of hybridization between the electron pockets, i.e |∆ in/out | > λ SOC . The first peak in the LDOS is due to the intra-band gaps ∆ in and ∆ out . We calculate δρ − (ω) and present the results in Fig.5(b). δρ − (ω) has positive peaks at ∼ 7 mev and ∼ 9 meV with no sign change between them indicating the sign change between ∆ in and ∆ out . Yet the third peak at ∼ 14 meV comes from inter-band pairing and causes two sign changes between second and third peak which is at odds with the experimental data 28 .

Spin-triplet dominated B2g-state
Even though a triplet dominated d-wave scenario is not realized in the phase diagram we quickly discuss this case with respect to criteria i)-iii) and QPI data.
As shown in Fig.6(a) and Fig.6(d) the LDOS is fully gaped and has a peak at 8 meV and 14 meV, moreover, the lower superconducting band shows back bending. ∆ in and ∆ out are plotted in Fig.6(c). Both gaps are equal in magnitude and nodeless due to smallness of hybridization between the electron pockets. In contrast to the singlet dominated case for the triplet dominated case we use |∆ B t | > |∆ B 1/3 | and sign(∆ B t ) = sign(∆ B 1 ) = sign(∆ B 2 ). In Fig.6(b) we plot δρ − (ω). Two small peaks between 7 meV and 8 meV reflect ∆ in and ∆ out contributions to δρ − (ω). The peak at 15 meV is due to the large inter-band gap. Overall, the behavior of δρ − (ω) is nearly consistent with experimental data as δρ − (ω) does not change sign in the experimentally relevant energy range. However, one has to bear in mind that in contrast to the triplet driven A 1g case the phase structure of the gaps on the Fermi surface affect the behavior of δρ − (ω) only near the lower peak, and the structure near the second (interband) peak is determined by the interband gap and does not necessarily bear information about the signs of the order parameters (see Appendix C.2). As we clearly see, to obtain the agreement with QPI experiments 28 the spin-triplet interband component from E g state is necessarily required in both symmetry states and triplet dominated A 1g state appears most likely one. This poses an important question how to detect it. To stay within QPI we propose to employ spin-resolved STM as an additional tool to further specify the underlying pairing state and also to distinguish between B 2g and A 1g -states.

E. Spin resolved STM
In the previous section we investigated corrections to the LDOS which arise from scattering on a single non-magnetic charge impurity and found the triplet dominated A 1g and B 2g pairing states in the presence of SOC to be consistent with the available experiments. As discussed in the introduction, however, the more conventional spin singlet A 1g and B 2g states have also been proposed. This raises the question whether SOC effects and the spin triplet order parameter can be further verified in experiment. Here we propose spin-resolved QPI as an additional tool to further specify the underlying pairing state. To illustrate the possible capabilities of this technique, we study the corrections to the LDOS in Born approximation from a single magnetic impurity which, within our basis of eq.(23), is given bŷ U = τ 3 ⊗ σ z ⊗V . Here τ and σ matrices act on Nambu and spin-space, respectively, and the matrixV describes orbital scattering. The Fourier transform of the spin resolved σ i -projected correction to the LDOS is given by [34] with the real and the imaginary part being even and odd in q, respectively. Note that a non-vanishing imaginary part requires the scattering potential to break inversion symmetry, which, as we show below can result from interorbital scattering. The scattering matrixV =V z +V x can be written aŝ where, in Born approximation,V z potentially contributes to the σ z -andV x to the σ x(y) -polarized LDOS when taking the trace in eq.(32). The labels intra and inter denote intra-and interorbital scattering, respectively. Note that in each block of the matricesV z (V x ) it is the off-diagonal (diagonal) elements that break inversion symmetry since under inversion Ψ X (−k) = σ z Ψ X (k) and Ψ Y (−k) = −σ z Ψ Y (k). Those elements contribute to Imδρ σ i (q, ω), while inversion symmetric ones contribute to Reδρ σ i (q, ω).
Results are presented in Fig.7 exemplary for ω = 8 meV where we show real and imaginary part of δρ σ i (q, ω) as a function of the scattering vector q for a σ z impurity with i = z (a-d) and i = x (e-h). Note that a similar behavior of these quantities is found for all frequencies within the range of the two maxima in the LDOS at 8meV and 14meV, respectively.' In Fig.7.(a,b) Reδρ σ z shows an even in q C 4 -symmetric QPI pattern qualitatively similar for the A 1g and B 2g pairing state. In contrast to that, see Fig.7.(c,d), Imδρ σ z is odd in q and exhibit a few qualitative differences between A 1g and B 2g especially along the intensity edges. Note, however, that Imδρ σ z is caused by the inter-orbital terms in V z and hence is expected to show a weakened intensity if V inter ≪ V intra .
In Fig.7.(e-h) we present real and imaginary part of the σ x -polarized QPI pattern δρ σ x (q, ω). We would like to stress that this quantity is non-zero only in presence of either a triplet gap ∆ A(B) t or SOC, which then itself induces a triplet gap. The same is true for the σ y -polarized state which is related to δρ σ x (q, ω) by 90 • rotation.
In Fig.7(e-f) we show the real part δρ σ x which is C 2 -symmetric and even in q while the imaginary part is odd, see Fig.7(g,h). In contrast to δρ σ z the patterns for A 1g and B 2g representation show qualitative differences in both real and imaginary part. Moreover, if the impurity occupies one Fe site (either on sub lattice A or B)V x contains scattering between xy X and xy Y components which breaks inversion symmetry and according to eq.(10) is of intraorbital nature. Hence, we expect the QPI signal of Imδρ σ x to be of the same order of magnitude as Reδρ σ z and hence may provide a quantity that allows to distinguish triplet driven A 1g from B 2g .
Thus we have shown that the spin-resolved QPI is a useful tool to study the details of the spin-orbit coupling and the triplet order parameter; our results suggest that a combination of spin-resolved QPI with the more conventional probes can be used to determine the symmetry as well as spin structure of pairing in iron-based superconductors.

V. CONCLUSION
In conclusion, we remind the reader that the highly electron-doped Fe-chalcogenides, in particular the (Li 1−x Fe x )OHFeSe system, have been discussed intensively in the context of the standard model of spin-fluctuation induced spin singlet pairing in the limit of weak spin-orbit coupling. In this case, the nodeless d-wave (B 2g ) and bonding-antibonding s ± (A 1g ) spin singlet states on the electron pockets have been considered most favorable, with the additional possibility of incipient s ± pairing via coupling to the incipient hole pocket at Γ. One of these states may eventually prove to be the correct pairing state for these materials, if the effective attraction in these channels generated by spin fluctuations driven by interpocket repulsion dominate the pair vertex 11 .
Here, we have instead studied the new possibility raised previously 20 , namely that interorbital triplet components generated by intrinsic attractions possible only in the presence of spin orbit coupling are responsible for the pairing in these systems. We explored the various possible Cooper-pairing symmetries using a mean-field decomposition of the Hubbard-Kanamori Hamiltonian including A 1g and B 2g states. For nonzero spin-orbit coupling, the superconducting order parameter is a combination of spin singlet and spin triplet gaps in each state. Treating attraction in singlet and triplet channels on equal footing and solving the selconsitency equations, we found for weak spin-orbit coupling a dominant spin singlet and small spin triplet gap yielding a state essentially equivalent to those identified in the usual spin fluctuation approach. For stronger spin-orbit coupling, however, the superconducting order parameter is a combination of spin singlet and dominant spin triplet gaps in each state. Focusing on the (Li 1−x Fe x )OHFeSe system, we identified the even parity A 1g -and B 2g -pairing states with a dominant spin triplet component to be consistent with available experiments, including current quasiparticle interference data, whereby according to our phase-diagram the A 1g state is slightly favored. The spin-singlet dominated A 1g and B 2g -states in this scenario without strong spin fluctuations are not consistent with at least one of the existing experiments.
In summary, to obtain a full moderately anisotropic gap on the Fermi pockets and its sign-changing character in agreement with experimental results on (Li 1−x Fe x )OHFeSe, we require either the traditional intraband A 1g or B 2g states generated by spin fluctuations, or the new triplet interband pair states in the same symmetry channels generated by intrinsic attraction in multiorbital correlated models. It is clearly of interest to identify experimental tests to distinguish between these possibilities. To this end we proposed using spin-polarized QPI to identify the possible triplet components present in the more exotic alternative states, and presented results for each of the triplet dominated states. A clear identification of triplet interband pairing using these results would be an important step forward in understanding the unusual superconductivity in the Fe chalcogenides. Blue curves refer to the same bands for zero superconducting gap. The low-lying energy band has its minima centered above the Fermi-level "back-bending".
In Fig.8(a) and Fig.8(d) we present the DOS and the energy band dispersion for the p-wave state. The DOS show two peaks tuned to lie at 8 meV and 14.3 meV, respectively and the lower energy band shows the back bending. Fig.  8(c) shows the absolute values of ∆ in and ∆ out projected on the Fermi surface pockets. In Fig.8 Here, one finds that the main feature of E u ⊗ U (1) in δρ − (ω) is that there are two sign changes between the peak energies of the total DOS (i.e. two negative peaks and positive values in between them), which does not agree with experiment 28 . Therefore this state does not seem consistent with the QPI data at least within Born scattering limit 28 . For consistency the spin-resolved QPI for the odd parity E u ⊗ U (1) p-state is presented in Fig.9. The SOC induced coupling between the spin singlet and triplet pairing channels translates into a coupling between intra-and inter-band order parameters in band space. It has been argued 21 that the second peak, seen in STM 28 can be at least partially due to an interband gap. In order to investigate how the predictions of HAEM's theory are affected by inter-band pairing we consider a simple model of a two band superconductor with superconductivity driven by intra-and inter-band pairing.
. The band dispersions are assumed to be simple parabolic ones ξ 1 (k) = k 2 2m − µ 1 and ξ 2 (k) = k 2 2m − µ 2 . Hence, ξ 2 = ξ 1 − 2B with B = 1 2 (µ 2 − µ 1 ) > 0. We only distinguish between intraand inter-band pairing and neglect for a moment their spin symmetry (i.e. consider them to be spin singlet) to get a better understanding on how the QPI data is affected by the inter-band pairing. We linearize both bands, and calculate δρ − (ω) and investigate how it depends on ∆ 12 , the relative phase between ∆ 12 and ∆ 21 and the band offset (direct gap).
For ω > 0 one needs to consider the cases: i) B < |∆ 12 | and ii) B > |∆ 12 | where in the former ω 1 and ω 3 and in the latter ω 2 and ω 4 correspond to the peak energies in δρ − (ω).

i) B < |∆|
In Fig.10, we plot δρ − (ω) in arbitrary energy units as a function of ω and B < |∆|. For ∆ 12 = ∆ 21 the behavior of δρ − is what we define here as "odd", meaning that between the two intraband gap energies δρ − changes sign, so the QPI-pattern is s ++ -like (solid blue curve). In case of ∆ 12 = −∆ 21 δρ − is "even", and one obtains as "s +− -like" pattern (solid orange curve). Hence, if the inter-band gap is larger than the band offset HAEM is sensitive to the relative phase between ∆ 12 and ∆ 21 .
ii) B > |∆| If B > |∆ 12 | the situation is a different one as can be seen by the dotted curves in Fig.10. In both cases ∆ 12 = ±∆ 21 and δρ − is "odd", mimicking an s ++ -pattern. Consequently, HAEM's method is not sensitive to the relative phase between ∆ 12 and ∆ 21 . Therefore, in a system with a large dominant inter-band gap the QPI signal depends not only on the relative phase between ∆ 12 and ∆ 21 but also on the ratio between |∆| and B, hence on the band structure.

Inter+Intra-Band
The influence of the inter-band gap on the QPI does affect the HAEM results on sign-changing and sign preserving intra-band gaps. To show this we numerically present the δρ − (ω) for ∆ 11 = ∆ 22 and sign(∆ 11 ) = −sign(∆ 22 ) leading to an initial s ± pattern. Then we increase |∆ 12 | and show that at a certain magnitude the phase information is lost. The results are plotted in Fig.11 where in the left and right panel we have ∆ 12 = −∆ 21 and ∆ 12 = ∆ 21 , respectively. Our results for the LDOS are presented in Fig.11(a) and Fig. 11(b). If ∆ 12 = 0, see orange curve, only intra-band If we switch the inter-band paring ∆ 12 = 1 on, see black curves, peak 3 and 4 appear as a consequence of inter-band pairing, whereas peak energies 1 and 2 are only slightly affected. If ∆ 12 is increased further, see blue curves, peak 2 and 3 merge.
In the absence of inter-band pairing the QPI pattern in Fig.11(c)-(d) is of s ± type, see orange curve. For moderate values ∆ 12 = 1 one can nicely distinguish between intra-and inter-band contributions Peak 1 and 2 show an s ± pattern due to the sign change between ∆ 1 and ∆ 2 , whereas between peak 3 and 4 the pattern is an s ++ since the inter-band gap is smaller that the band offset B, see Fig.10. For moderate values of ∆ 12 , see blue curves, peak energies 1 and 2 start to deviate from ∆ 1 and ∆ 1 as marked by the vertical dashed lines. As soon as peak 2 and 3 merge the s ± pattern between peak 1 and peak 2 becomes more and more difficult to resolve. The blue curve corresponds to the situation where the inter-band gap ∆ 12 is larger than the intra-band gap but smaller that the band offset. The latter condition causes an s ++ pattern between peak 3 and 4 independent of the relative phase between ∆ 12 and ∆ 21 , sec.C 1. This is accompanied by ∆ 11 and ∆ 22 having a s ± symmetry. The combination of both patterns leads to a sign-change of δρ − (ω) in the region 2 < ω < 4 which might be misinterpreted as a an s ++ pattern between ∆ 11 and ∆ 22 and hence carries no concrete phase information.